{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 机器学习练习 - 逻辑回归"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "代码更新地址：https://github.com/fengdu78/WZU-machine-learning-course\n",
    "\n",
    "代码修改并注释：黄海广，haiguang2000@wzu.edu.cn\n",
    "\n",
    "在这一次练习中，我们将要实现逻辑回归并且应用到一个分类任务。我们还将通过将正则化加入训练算法，来提高算法的鲁棒性，并用更复杂的情形来测试它。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 逻辑回归"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "在训练的初始阶段，我们将要构建一个逻辑回归模型来预测，某个学生是否被大学录取。设想你是大学相关部分的管理者，想通过申请学生两次测试的评分，来决定他们是否被录取。现在你拥有之前申请学生的可以用于训练逻辑回归的训练样本集。对于每一个训练样本，你有他们两次测试的评分和最后是被录取的结果。为了完成这个预测任务，我们准备构建一个可以基于两次测试评分来评估录取可能性的分类模型。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "让我们从检查数据开始。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 29,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import pandas as pd\n",
    "import matplotlib.pyplot as plt"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 30,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>Exam 1</th>\n",
       "      <th>Exam 2</th>\n",
       "      <th>Admitted</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>34.623660</td>\n",
       "      <td>78.024693</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>30.286711</td>\n",
       "      <td>43.894998</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>35.847409</td>\n",
       "      <td>72.902198</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>60.182599</td>\n",
       "      <td>86.308552</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>79.032736</td>\n",
       "      <td>75.344376</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "      Exam 1     Exam 2  Admitted\n",
       "0  34.623660  78.024693         0\n",
       "1  30.286711  43.894998         0\n",
       "2  35.847409  72.902198         0\n",
       "3  60.182599  86.308552         1\n",
       "4  79.032736  75.344376         1"
      ]
     },
     "execution_count": 30,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "path = 'ex2data1.txt'\n",
    "data = pd.read_csv(path, header=None, names=['Exam 1', 'Exam 2', 'Admitted'])\n",
    "data.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 31,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(100, 3)"
      ]
     },
     "execution_count": 31,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "data.shape"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "让我们创建两个分数的散点图，并使用颜色编码来可视化，如果样本是正的（被接纳）或负的（未被接纳）。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 32,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAtMAAAHgCAYAAABn8uGvAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAgAElEQVR4nOzdf3hcZZ3//9c7FCgmowhWv1iEIo0/kB8F4w80HyJUUQtIWtGgfhXdrnX94Nofutu63+tS18uPKPrdGr7uR5ctC7jrlgibBK616y+UDxZddlNEQcBNdQvWoMUKMgmiQN7fP+45zSSdSSYnM3POmfN8XFeuk7nPTObOyWTmfe7zvt+3ubsAAAAAzF9b0h0AAAAAsopgGgAAAIiJYBoAAACIiWAaAAAAiIlgGgAAAIiJYBoAAACIaVHSHViIZz3rWb5s2bKkuwEAAIAWt2vXrt+4+5KZ7ZkOppctW6aRkZGkuwEAAIAWZ2b3V2onzQMAAACIiWAaAAAAiIlgGgAAAIgp0znTAAAAWfXEE09o7969evzxx5PuCsosXrxYxx57rA499NCa7k8wDQAAkIC9e/eqUCho2bJlMrOkuwNJ7q79+/dr7969OuGEE2p6DGkeAAAACXj88cd19NFHE0iniJnp6KOPntfVAoJpAACAhBBIp898/yYE0wAAADk2NDQkM9N9991Xcf+73/1u3XDDDTX/vLGxMV100UWSpDvvvFM7duw4sO+WW27R97///Xn3cdmyZfrNb34z78c1A8E0AABABhSL0rZt0ubNYVss1ufnbt++Xd3d3bruuuvq8vOe+9znHgi+6xVMpxnBNAAAQMrt3CktXSpt2CBdfnnYLl0a2hdifHxct912m6666qoDwbS76wMf+IBOOukknXfeedq3b9+B+y9btkx/9Vd/pTPPPFNdXV2644479PrXv14nnniivvSlL0mS9uzZo5NPPll//OMf9dGPflQDAwNasWKFPvOZz+hLX/qStm7dqhUrVuh73/ueHnroIb35zW/Wy172Mr3sZS/TbbfdJknav3+/zj33XJ1++ul63/veJ3df2C/aQA2r5mFm/yDpfEn73P3kUttRkgYkLZO0R9Jb3f1hC8kp/ZJWSXpM0rvd/Y5G9Q0AACArikVp1arpI9ETE2G7apU0NiZ1dMT72cPDw3rDG96gF7zgBTrqqKN0xx13aM+ePfrpT3+qu+66S7/+9a910kkn6U/+5E8OPOZ5z3uefvCDH2jjxo1697vfrdtuu02PP/64XvKSl+jP/uzPDtzvsMMO0yc+8QmNjIzoC1/4giTp97//vTo6OvThD39YkvT2t79dGzduVHd3tx544AG9/vWv17333qu//uu/Vnd3tz760Y/qa1/7mq688sp4v2ATNLI03jWSviDpy2VtWyTd7O6fNrMtpdubJb1RUmfp6xWSvljaAgAA5NrAgDQ5WXnf5GTYv3ZtvJ+9fft2bdiwQZJ08cUXa/v27XriiSf0tre9TYcccoie+9zn6pxzzpn2mDe96U2SpFNOOUXj4+MqFAoqFApavHixHnnkkXk9/7e//W3dc889B24/+uijKhaLuvXWWzU4OChJOu+88/TMZz4z3i/YBA0Lpt39VjNbNqP5QkmvKX1/raRbFILpCyV92cMY/r+b2ZFmdoy7P9io/gEAAGTB6OjUSPRMExPS7t3xfu7+/fv1ne98R3fffbfMTE899ZTMTKtXr561osXhhx8uSWprazvwfXT7ySefnFcfJicn9YMf/EBHHHHEQfuyUumk2TnTz4kC5NL22aX2pZJ+UXa/vaU2AACAXOvslNrbK+9rb5eWL4/3c2+44Qa9613v0v333689e/boF7/4hU444QQdddRRuu666/TUU0/pwQcf1He/+93YfS8UCiqW5afMvH3uueceSAGRwoRFSTrrrLP0la98RZL0b//2b3r44Ydj96HR0jIBsdKpR8VMczNbZ2YjZjby0EMPNbhbAAAAyerrk9qqRGxtbWF/HNu3b9fq1auntb35zW/Wr371K3V2duqUU07R+9//fvX09MR7Aklnn3227rnnHq1YsUIDAwO64IILNDQ0dGAC4hVXXKGRkRGdeuqpOumkkw5MYvzYxz6mW2+9VWeccYa++c1v6rjjjovdh0azRs6OLKV5/GvZBMSfSnqNuz9oZsdIusXdX2hmf1f6fvvM+83287u6unxkZKRh/QcAAGiUe++9Vy9+8Ytruu/OnWGy4eRkSO1obw+B9I4dUnd3gzuaQ5X+Nma2y927Zt632SPTN0m6pPT9JZJuLGt/lwWvlPQ78qXn4C4NDYVtLe0AACCzurtD1Y7+fmnLlrAdGyOQToNGlsbbrjDZ8FlmtlfSxyR9WtJXzWytpAckvaV09x0KZfF2K5TGe0+j+tUyhoelNWuk9eulrVslsxBAb9wY/sMGB6UZl24AAEB2dXTEr9qBxmlkNY+3Vdm1ssJ9XdKljepLS+rtDYF0f3+4vXXrVCC9fn3YDwAAgIZqZJ1pNJJZCKClEEBHQXX5SDUAAAAaKi3VPBBHeUAdIZAGAABoGoLpLItypMtt3MjkQwAAgCYhmM6q8smG69eHWjlRDjUBNQAAqIGZ6UMf+tCB25/73Of08Y9/fNbHDA8PT1sCvJLTTjtNb3tbtelz0p49e3TyySfPq68f/ehH9e1vf1uS9PnPf16PPfbYgX2f+tSn5vWzJOmaa67RBz7wgXk/biaC6awaHp4KpKPUjq1bpwLq4eGkewgA81IsStu2SZs3h23ZImlAvjWwHO7hhx+uwcFB/eY3v6n5MXMF0/fee68mJyd16623aqLaOugxfOITn9BrX/taSfUJpuuFYDqrentD+bvyHOkooB4cbKlqHnzAAq1v505p6VJpwwbp8svDdunS0A7kXlQOt/zKc3SFes2aBQ2gLVq0SOvWrdPWmXOwJN1///1auXKlTj31VK1cuVIPPPCAvv/97+umm27SX/zFX2jFihX62c9+dtDj/vmf/1nvfOc7de655+qmm2460L5r1y6ddtppOvPMM/W3f/u3B9qvueYa9fb26oILLtAJJ5ygL3zhC/qbv/kbnX766XrlK1+p3/72t5Kkd7/73brhhht0xRVXaGxsTGeffbbOPvtsbdmyRb///e+1YsUKveMd75Ak/dM//ZNe/vKXa8WKFXrf+96np556SpJ09dVX6wUveIF6enp02223xT5u07h7Zr9e+tKXOlrb977nXii4t7e7S2FbKIR2AK3h0UfD/3WIDqZ/FQruxWLSPQQa45577qntjpOT7uvXh3+K9esr346pvb3df/e73/nxxx/vjzzyiH/2s5/1j33sY+7ufv755/s111zj7u5XXXWVX3jhhe7ufskll/j1119f9Wd2dnb6nj17/Bvf+IZfcMEFB9pPOeUUv+WWW9zd/cMf/rC/5CUvcXf3q6++2k888UR/9NFHfd++ff70pz/dv/jFL7q7+4YNG3zr1q0HPe/xxx/vDz300LTfI3LPPff4+eef73/84x/d3f3973+/X3vttT42NubPe97zfN++ff6HP/zBX/WqV/mll15a8Xeo9LeRNOIV4lFGppFaxWJYOrVYDEunSmEbtY+PJ9s/APUxMBCmfVQyORn2A7k2M5Wzre3gVM8FePrTn653vetduuKKK6a1/+AHP9Db3/52SdI73/lO7azhUtF//ud/asmSJTr++OO1cuVK3XHHHXr44Yf1u9/9To888oh6enoO/LxyZ599tgqFgpYsWaJnPOMZuuCCCyRJp5xyivbs2TOv3+fmm2/Wrl279LKXvUwrVqzQzTffrJ///Oe6/fbb9ZrXvEZLlizRYYcdpr6+vnn93GoIppFafMAC+TA6OnXCPNPEhLR7d3P7A6RSg8vhbtiwQVddddWsOc5Ww3Nt375d9913n5YtW6YTTzxRjz76qP7lX/5F7j7r4w8//PAD37e1tR243dbWpieffHIev0nIurjkkkt055136s4779RPf/rTA5Mqa/kd5otgGqnFByzyJq/zAzo7pfb2yvva26Xly5vbHyCVGlwO96ijjtJb3/pWXXXVVQfaXvWqV+m6666TJH3lK19Rd3e3JKlQKKhY4Q1qcnJS119/vX784x9rz5492rNnj2688UZt375dRx55pJ7xjGccGN3+yle+sqD+zuzDoYceqieeeEKStHLlSt1www3at2+fJOm3v/2t7r//fr3iFa/QLbfcov379+uJJ57Q9ddfv6A+RAimkVp8wCJP8jwBr68vXLWupK0t7Mfc8noylgtNKof7oQ99aFpVjyuuuEJXX321Tj31VP3jP/6j+kurLV988cX67Gc/q9NPP33aBMRbb71VS5cu1dKlSw+0nXXWWbrnnnv04IMP6uqrr9all16qM888U0ccccSC+rpu3Tq98Y1v1Nlnn33g9qmnnqp3vOMdOumkk/TJT35S5557rk499VS97nWv04MPPqhjjjlGH//4x3XmmWfqta99rc4444wF9SFinuF6xF1dXT4yMpJ0N9AgxWIIJip9IBQK0tiY1NHR/H4BcRWLIT1pdDScLPb1hdcyr/Vw0rBqVYgRJibCCXNbm7Rjh1QaDMMsOH7ZdO+99+rFL37x3HccGgpVO8pzpMsD7MFBafXqxnc4Ryr9bcxsl7t3zbzvoqb1CpinQiF8EFT7gGj14AKtpVKws2lTeC3fd9/c8wPWrm1uf5utuzucNAwMhBSu5cvDyQb/53Mrn6wdiVLkVq3Kx8lYy4vK4fb2HlwOt6enpcrhZhHBNFKND1i0grmCnT/9U+YHSOH/uhknDdWuEGRVLZO1W/1krOWZVR55rtaOpiKYRuo16wMWaJS5gp39+8NIdaWAmvkB9TXbFYKspkMwWRtIFhMQAaDB5gp2jj6aCXjN0Kq165msnW1ZnrvWqub7NyGYTgv3MMFg5h+wWjuAzJgr2HnJS8LIaKEwdb/29ql5A6Q11Uer1q6nGkp2LV68WPv37yegThF31/79+7V48eKaH0OaR1oMDzNTF2hRfX0hlaCSKNjp6GB+QKO1ajoEk7Wz69hjj9XevXv10EMPJd0VlFm8eLGOPfbYmu9PMJ0Wvb1TNSOlEFCX15Rkpi6QWbUGO8wPaKzoCkEr5qYzWTubDj30UJ1wwglJdwMLRJ3pNCkfiY6Uj1QDyLTxcYKdJFHPG8BCVKszTTCdNu7Tk98mJwmkAaBOWNwEQFws2pIF0ch0uY0bGZkGgDohHQJAvRFMp0V5ikeU2lGe8kFADQB1QW46gHoimE6L4eHpgXS0TKgU2nt6qOYBAACQMgTTadHbG8rf9fZOjUBHAXVPD9U8AAAAUohgOi3MKo88V2sHAABA4lgBEQAAAIiJYBoAAACIiWAaAAAAiIlgGgAAAIiJYBoAAACIiWoeAIDMKRbDKoajo1JnZ1jFsFBIulcA8ohgGgCQKTt3SqtWSZOT0sSE1N4ubdok7dgRlgsHgGYizQMAkBnFYgiki8UQSEthG7WPjyfbPwD5QzANAMiMgYEwIl3J5GTYDwDNRDANAMiM0dGpEemZJiak3bub2x8AIGcaAJAZnZ0hR7pSQN3eLi1f3vw+ZR2TOYGFMXdPug+xdXV1+cjISNLdAAA0SbEoLV0atjMVCtLYmNTR0fx+ZVWlyZxtbUzmBCoxs13u3jWznTQPAEBmFAoh0CsUQuAnhW3UTiBdOyZzAvVBmgcAIFO6u8MI9MBAyJFevjykJhBIz08tkznXrm1un4AsIpgGAGRORweB3kIxmROoD9I8AADIoWgyZyVM5gRqRzANAEAO9fWFyYaVtLWF/QDmRjANAEAOlU/mfNrTQtuiRdLhh0s33EAOOlArgmkAAHKquzsEzpOT0qGHSk8+GQLqiy4KZfMAzI1gGgCAnCoWQ+D8+OPSE0+ENsrjAfNDMA0AQIMUi9K2bdLmzWFbabGZJNVSHg/A7CiNBwBAA1RaXXDTpnStLkh5PGDhEhmZNrP1Zna3mf3EzDaU2o4ys2+Z2Whp+8wk+gYAwEJlZXVByuMBC9f0YNrMTpb0Xkkvl3SapPPNrFPSFkk3u3unpJtLtwEAyJyspE9QHg9YuCRGpl8s6d/d/TF3f1LS/5G0WtKFkq4t3edaSb0J9A0AgAXLSvpEeXm8aIS6vX2qnfJ4wNySyJm+W9L/MrOjJf1e0ipJI5Ke4+4PSpK7P2hmz06gbwAAzFuxGEabR0dD6sTznheC0koBddrSJ7q7pbGx0P/du0Pf+voIpIFaNT2Ydvd7zewzkr4laVzSjyQ9WevjzWydpHWSdNxxxzWkjwAA1KrSREOz6mkeaUyf6OiQ1q5NuhdANiUyAdHdr3L3M9z9LEm/lTQq6ddmdowklbb7qjz2SnfvcveuJUuWNK/TAOom7eXCgFpVm2g4Ph4C6o4O0ieAVpdIaTwze7a77zOz4yStkXSmpBMkXSLp06XtjUn0DUBjZaFcGFCr2SYamkmf/rS0eDHpE0ArS6rO9L+UcqafkHSpuz9sZp+W9FUzWyvpAUlvSahvABqkfBQvEo3mrVoV8jYJNJAlc0003LtXuuyy5vYJQHMlEky7+/+o0LZf0soEugOgSWopF0beJrIkqtOchYmGSLeZk1j7+kJaENKPFRCBFGn1N9OslAsDatXXF9KUKknjREOkE+lv2UYwDaREHt5MGcVDq4kmFM78321rY6IhakP6W/YlUs0DwHRZWXp4oVhtDa0oqtPc3y9t2RK2Y2OtcxKMxsrKapmojpFpIAXykkvMKB5aFXWaERfpb9lHMA2kQJ7eTFltDQCmkP6WfQTTQArk7c2UUTwAeTLb5HImsWafuXvSfYitq6vLR0ZGku4GsGDForR0aeWVAAsFJqAAraTVq/ZgukqTy6PUtiivvpb7IHlmtsvduw5qJ5gG0oE3U6D18X+eL/MZKBkfJ/0t7aoF06R5oKVlaQSIXGKgtVECLX/mM7mc9LfsIphGy8pi3WbeTIHWlZeqPZiSp8nleUYwjZbECBCAtCGwyp+8TS6vJEtXiONi0Ra0JIrgA0ibKLCqJC+BVd7kfaGqnTtDzviGDdLll4ft0qWhvZUQTKMlMQKEJBSL0rZt0ubNYVtp0hHyK++BVR5FC1UVClMnUu3tU+2tfIU0Lyv7SqR5oEVxaQ3NlsUcfTQXK4DmU14nl+dpjgDBNFoSRfDRTOToo1Z5DazyLo+Ty/N0hZhgGi2JESA0U55GYFC7ahOv8hhYIX/ydIWYYBotixEgNEueRmAWIg+z+iOk/SDv8nSFmGAaLY0RIDRDnkZg4spTcEnaD5CvK8RU8wCABaJKw+zyNKtfojQnEImuEPf3S1u2hO3YWOudQDMyDQALlKcRmDjyllNO2g8wJQ9XiAmmAaAOyNGvLm/BJWk/QL4QTANAneRhBCaOVg0uq02ozNPEKwCSuXvSfYitq6vLR0ZGku4GAGAWxWJYQrjSipCFQjYn5FWaUBml9XR3z70fQPaY2S5375rZzsg0AKChWi2nvJZqHaT9APlBMA0AaLhWCi5rnVBJ2g+QDwTTAICmaJXgMm8TKgHMjjrTAADMQzShspIsT6gEEA/BNAAA88AiPQDKEUwDADAP0YTKQmFqhLq9fao9i3ngAOIjZxoAgHlqpQmVABaGYBoAgBhaZUIlgIUhzQMAAACIiWAaAAAAiIlgGgAAAIiJYBoAAACIiWAaAAAAiIlgGgAAAIiJYBoAAACIiWAaAAAAiIlgGgAAAIiJYBoAAACIiWAaAAAAiIlgGgAAAIhpUdIdAADUX7EoDQxIo6NSZ6fU1ycVCkn3CgBaD8E0ALSYnTulVaukyUlpYkJqb5c2bZJ27JC6u5PuHQC0FtI8AKCFFIshkC4WQyAthW3UPj6ebP8AoNUQTANACxkYCCPSlUxOhv0AgPohmAaAFjI6OjUiPdPEhLR7d3P7AwCtLpFg2sw2mtlPzOxuM9tuZovN7AQzu93MRs1swMwOS6JvAJBlnZ0hR7qS9nZp+fLm9gcAWl3Tg2kzWyrpg5K63P1kSYdIuljSZyRtdfdOSQ9LWtvsvgFA1vX1SW1V3tnb2sJ+AED9JJXmsUjSEWa2SNLTJD0o6RxJN5T2XyupN6G+AUBmFQqhakehMDVC3d4+1d7RkWz/AKDVNL00nrv/0sw+J+kBSb+X9E1JuyQ94u5Plu62V9LSZvcNAFpBd7c0NhYmG+7eHVI7+voIpAGgEZoeTJvZMyVdKOkESY9Iul7SGyvc1as8fp2kdZJ03HHHNaiXAJBtHR3SWpLlAKDhkkjzeK2k/3b3h9z9CUmDkl4l6chS2ockHStprNKD3f1Kd+9y964lS5Y0p8cAAABABUkE0w9IeqWZPc3MTNJKSfdI+q6ki0r3uUTSjQn0DQAAAKhZ04Npd79dYaLhHZLuKvXhSkmbJW0ys92SjpZ0VbP7hhbhLg0NhW0t7QAAADElUs3D3T/m7i9y95Pd/Z3u/gd3/7m7v9zdl7v7W9z9D0n0DS1geFhas0bauHEqcHYPt9esCfsBAADqoOkTEIGG6+2V1q+X+vvD7a1bQyDd3x/ae6m6CORFsRiqmoyOhgVt+vpCmUAAqBfzDF/y7urq8pGRkaS7gTSKRqKjgFoKgfTWrZJZcv0C0DQ7d0qrVkmTk2Ep9fb2sHDNjh2hfCAAzIeZ7XL3roPaCaZbjHtIY+jtnR40VmtvZe7Tl4KbnMzP745cYzQ2HIOlS8N2pkIh1OGm7jaA+agWTCe1AiIahXzhIPqdy5UfE6BF7dwZgsgNG6TLLw/bpUtDe54MDITz50omJ8N+AKgHgulWU54vHAWPecsXnvk7T04efEyAFlQshrSGYjGkNUhhG7WPjyfbv2YaHZ06BjNNTISVIQGgHpiA2GrMQl6wFILHKGc4T/nCw8NTgXT0O5cfk54eafXqZPsINEAto7F5WRWxszPkSFcKqNvbwxLrAFAPjEy3ovLgMZKXQFoKo++Dg9N/5+iYDA7mY3QeucRo7JS+vulTJsq1tYX9AFAPBNOtKO/5wmZh5HnmyUO1dqBFRKOxleRtNLZQCFU7CoWpY9LePtXO5EPkRbEobdsmbd4ctpUm5WJhqObRambmC8+ssZynEWogZ6hgcbDx8ZDesnt3OJno68vfMUB+UR6yviiNlxdDQ6FqR3ngXB5gDw6SLwy0MD48AUicXDdCtWCaCYitJsoXLq8nHeUL9/SQLwy0uO7u8CHJaCyQb0xIbh6C6VYT5QXX2g6g5XR08CGJ5mOxoHRhQnLzEEwDAIAFqZRetGkT6UVJojxk81DNAwAAxMZiQelEecjmIZhGctzDhMmZk2CrtQMAUoel29OJ8pDNQ5oHkjM8TOWRhJHjODuODzA3cnPTiwnJzUEwjeT09oZAOlryfGZNbCqPNBQ5jrPj+AC1ITc33ZiQ3HjUmUayykeiIywu03DUH50dxweoHf8vyItqdabJmc66rOcdRzWwyxFINxw5jrPj+KAVNGsZaXJzkXcE01kX5R1v3DgVOEejvWvWhP1pFvW1XPnvgoYgx3F2HB9k3c6dYbR4wwbp8svDdunS0N4IUW5uf7+0ZUvYjo2REoV8IGc667Kcd1ye4hGldpSnfDBC3TDkOM6O44MsKy9VF4ley6tWNS7tgtxc5BUj01kXpUlEAXVb2/TgNM3B6PDwwX0t/13SPqqeYdQfnR3HB1lGmhLQXATTrSCrece9vaH8XXlfo99lcDDdo+oZR47j7Dg+yDLSlIDmIs2jFVTLO057QG1WuY50tXbUFfVHZ8fxQVaRpgQ0F6Xxsm62vOMspHoAAOqKUnVAY1QrjcfIdNZVyzuWQntPD6O8AJAjUTrSzEWH2tqmpymxwidQH4xMZ517CKh7e6ePQFdrBwDkwvh49TSlSit8RsE25eyAyqqNTBNMAwCQI6SBAPGwAiIAAKB0HlBnBNMAAOQIpfOA+iKYBgAgR6LSeZVQOg+YP4JpAKizYlHatk3avDlsK+WmAklhhU+gvgim88xdGhoK21raAcxp584wuWvDBunyy8N26dLQHiHYRpJY4ROoL6p55NnQkLRmzfQa1eWLwAwOUqMamIdaqiTceSclyZAOs5XOA3AwFm3BwXp7QyDd3x9uz1w9sbc32f4BGTNXlYRrr5U+8pHpwXY0EWzVKkqSobk6OqS1a5PuBZB9BNN5NnO1xCioZhlyIJa5qiT867/OXZKM4AYADpbmFTvJmc678oA6QiANxDJXlQSJkmQAMF+1zEVJEsF03kU50uU2bmTyIRDDXFUSzjuPkmQAMB/FYkiDKxanBiMmJqbax8eT7Z9EMJ1v5ZMN168P15mjHGoCamDe5qqScMkllCQDgPnIwoqd5Ezn2fDwVCAdpXaU51D39FDNA5in7u4wkbBalYQdO6pX82DyIQBMl4UVOwmm86y3N5S/6+2dypGOAuqeHqp5ADHNViVhrmAbADAlmotSKaBOS3ocdaYBAACQSrXU72/WYES1OtPkTCN9WJkRAAAoGyt2EkwjfYaHw8qM5ZMgo8mSa9aE/QAAIBei9Lj+fmnLlrAdG0vPqrHkTCN9WJkRAACUSfOKnQTTSB9WZgQAABnBBESkl/v0oryTkwTSAHIjzcsnA3nEBERkCyszAsixtC+fDGAKwTTSh5UZAeRYFpZPBjCl6TnTZvZCSeWLPz5f0kclfbnUvkzSHklvdfeHm90/pAArMyLnuLyfbo3++9SyfHJaJ2K1Ov43UcmcOdNm9gJJX5T0HHc/2cxOlfQmd//kgp/c7BBJv5T0CkmXSvqtu3/azLZIeqa7b57t8eRMtyj3EFCXr8w4WzvQQnburL7ceFrKQOVZM/4+mzeH1I5qtmyRLrusPs+F2vG/iYXkTP+9pI9IekKS3P3Hki6uU79WSvqZu98v6UJJ15bar5VE/bO8MgsjzzMD5mrtQIvg8n66NevvEy2fXElalk/OG/43MZtagumnuft/zGh7sk7Pf7Gk7aXvn+PuD0pSafvsSg8ws3VmNmJmIw899FCdugGgVsWitG1bGD3btq3yEq+Ip5bL+0hOs/4+fX3TCxmVa2sL+9Fc/G9iNrUE078xsxMluSSZ2UWSHlzoE5piv4IAACAASURBVJvZYZLeJOn6+TzO3a909y5371qyZMlCuwFgHqgw0Fijo1OjXjNNTEi7dze3P5iuWX+fLCyfnDf8b2I2tUxAvFTSlZJeZGa/lPTfkt5Rh+d+o6Q73P3Xpdu/NrNj3P1BMztG0r46PAeAOim/zBmJPlxWrQpLu/IhvzDR5f1KH9pc3k9eM/8+0fLJAwMhUFu+PIxId3QwCS4J/G9iNrNOQDSzNkkXuftXzaxdUpu71+WirpldJ+kb7n516fZnJe0vm4B4lLv/5Ww/gwmIQPNs2xZGoqt9mPT3U2FgoYrFMNJfKXWmUOCEJWlp+PswCS4ZafjbZ9FCTvzSeNJYbQLirCPT7j5pZh+Q9FV3r3KBI1ZnnibpdZLeV9b8aUlfNbO1kh6Q9JZ6PR+AheMyZ+NFl/GrBUt8WCcr6b8PV4eSk/TfPosqnfht2lTbid9CHpuEWtI8vmVmH1aoAX3go9Tdfxv3Sd39MUlHz2jbr1DdA0AKcZmzOWa7vI/kJfn3of50svjfrN1CTvyyeNJYSzD9J6XtpWVtrrDYCoCc6OsLIwOVUGGgvjo6CIrSLKm/D1eHksf/Zm0WcuKXxZPGOYNpdz+hGR0BkG5c5gSSxdUhZMVCTvyyeNI4ZzBtZodKer+ks0pNt0j6O3d/ooH9ApBCXOYEksPVIWTFQk78snjSWMty4tskHaqp1QnfKekpd//TBvdtTlTzAADkCdU8kAULqX6S5sopsap5lLzM3U8ru/0dM/tR/boGAABqwdUhZMFC0gKzmFJYSzD9lJmd6O4/kyQze76kpxrbLQAAUAmT4JAFCznxy9pJYy3B9F9I+q6Z/VySSTpe0nsa2isAAABk2kJO/LJ00lhLNY+bzaxT0gsVgun73P0PDe8ZAAAAkHJtc93BzC6VdIS7/9jdfyTpaWb2PxvfNQAAACDd5gymJb3X3R+Jbrj7w5Le27guAQAAANlQS850m5mZl2romdkhkg5rbLcAAFlRLIaJQqOjoUZsX1+YkQ8AeVDLyPQ3JH3VzFaa2TmStkv6emO7BaSQuzQ0FLa1tAM5sHNnqAm7YYN0+eVhu3RpaAeAPKglmN4s6WaFVRAvLX3/l43sFJBKw8PSmjXSxo1TgbN7uL1mTdgP5EixGGrBFotTq5VNTEy1j48n2z8AaIY5g2l3n3T3L0l6u6RPShpyd+pMI396e6X166X+/qmAeuPGcHv9+rAfyJGBgbCoQiWTk2E/ALS6qjnTZvYlSf+fu//EzJ4h6QcKi7UcZWYfdvftzeokkApm0tat4fv+/vAlhUB669awH8iR0dGpEemZJibCYgtAXjB3IL9mG5n+H+7+k9L375H0X+5+iqSXijQP5FV5QB0hkEZOdXaGZX4raW8Pq5YBecDcgXybLZj+Y9n3r5M0LEnu/quG9ghIsyi1o1x5DjWQI319UluVT5G2trAfaHVZmjtQLErbtkmbN4dtsZh0j1rDbMH0I2Z2vpmdLunVKlXwMLNFko5oRueAVJmZIz05eXAONZAjhYK0Y0fYRiPU7e1T7R0dyfYPaIaszB1g9LxxZqsz/T5JV0j6vyRtKBuRXinpa43uGJA6w8NTgXSU2lGeQ93TI61enWwfgSbr7pbGxkLAsHt3SO3o6yOQRn5kYe5A+eh5JOrzqlXhf5j/2fiqBtPu/l+S3lCh/RsKtaeBubmHILS3d3pecbX2NOvtlQYHp/c5Cqh7eqjm0UBM7Em3jg5p7dqkewEkI5o7UCmgTsvcgVpGz/kfjq+WOtNAfK1Um9ksjDzPDP6rtaMuuDQJIM2yMHcgC6PnWUYwjcaiNjMWIEsTewDkUxbmDlB5p7Fmy5kGFo7azFgALk0CyIK0zx3o65M2baq8Ly2j51k2azBtZi+StFTS7e4+Xtb+Bnf/eqM7hxYRBdRRIC0RSKMmXJoEkBVpnjsQjZKvWhUGIiYmwoh0W1t6Rs+zrGqah5l9UNKNkv5c0t1mdmHZ7k81umNoIdRmRkxcmgSA+ohGz/v7pS1bwnZsLLRjYWYbmX6vpJe6+7iZLZN0g5ktc/d+SQwpojYzc6S3bp26LTFCjVlxaRIA6ifNo+dZNlswfUiU2uHue8zsNQoB9fEimEatqM2MBeDSJAAg7cyrXGo3s+9I2uTud5a1LZL0D5Le4e6HNKeL1XV1dfnIyEjS3cBsWqnONBIzPp7eiT0AgHwws13u3nVQ+yzB9LGSnixb+bB836vd/bb6d3N+CKYBAADQDNWC6dlWQNw7y77EA2kAAAAgaSzaAqB27tLQ0MGVWKq1AwDQ4gimAdSulZaHBwCgDmpeAdHMnl5+f3f/bUN6BCC9ypeHl6aXOmR5eABADs0ZTJvZ+yR9QtLvJUXXcF3S8xvYLwBpxPLwAABMU7Wax4E7mI1KOtPdf9OcLtWOah5AQtxDsefI5CSBNACgpVWr5lFLzvTPJD1W/y4ByCSWhwcA4IBacqY/Iun7Zna7pD9Eje7+wYb1CkA6sTw8AADT1BJM/52k70i6S9JkY7sDINVYHh4AgGlqCaafdPdNDe8JgPTr7ZUGB6cvAx8F1D09VPMAAOROLTnT3zWzdWZ2jJkdFX01vGcA0scsjDzPTOWo1g4AQIurZWT67aXtR8raKI0HAACA3JszmHb3E5rREQAAACBraloB0cxOlnSSpMVRm7t/uVGdAgCkU7EoDQxIo6NSZ6fU1ycVCkn3CkArydr7TC2LtnxM0msUgukdkt4oaae7X9Tw3s2BRVuAhLiHyh7lExFna0dL2LlTWrUqrNEzMSG1t4e1e3bskLq7k+4dgFaQ5veZhSzacpGklZJ+5e7vkXSapMPr3D8AWTI8LK1ZM32xlqgG9Zo1YT9aSrEYPuCKxfABJ4Vt1D4+nmz/AGRfVt9nagmmf+/uk5KeNLOnS9onJh+mk7s0NHTwSnTV2oG4entDren+/qmAunwxF0rktZyBgTBSVMnkZNgPAAuR1feZWoLpETM7UtLfS9ol6Q5J/9HQXiEeRgvRLFFt6Sigbms7eDEXtJTR0amRopkmJqTdu5vbHwCtJ6vvM3MG0+7+P939EXf/kqTXSbqklO6BtGG0EM1UvvphhEC6ZXV2htzFStrbpeXLm9sfAK0nq+8zcwbTZrY2+t7d90j6SWlSYmxmdqSZ3WBm95nZvWZ2ZmkxmG+Z2Whp+8yFPEcuMVqISDNSfqKTtXLlV0XQUvr6wltKJW1tYT8ALERW32dqSfNYaWY7Sisgnizp3yUttEBJv6Svu/uLFCY03itpi6Sb3b1T0s2l25gvRgshNT7lZ+ZVj8nJg6+KoKUUCmE2faEwNXLU3j7V3tGRbP8AZF9W32dqWbTl7WbWJ+kuSY9Jepu73xb3CUuTGM+S9O7Sz/+jpD+a2YUKJfgk6VpJt0jaHPd5cqvaaCEBdb6Up/xI4e9fz5Sf4eGDr3pEJ3H9/VJPT1hePCWyVrM0rbq7pbGxcCx37w6XXPv60vsBByB7svg+U0ud6U6F4PYuSS+WdI+kTe7+WKwnNFsh6crSzzlNYVLjekm/dPcjy+73sLsflOphZuskrZOk44477qX3339/nG60ppmjhTMDKALqfCl/PUTq9TrIUJ3pNNcsBQBkR7U607UE0/dJutTdbzYzk7RJ0p+4+0tidqRLIVXk1e5+u5n1S3pU0p/XEkyXY9GWGYaGwiX88oCpPKAaHEzVaCGawH16AtrkZGqC3GYoFqWlS8N2pkIhjH6kebQDAJAeC1m05eXufrMkefD/SlrINeK9kva6++2l2zdIOkPSr83smFJnj1GoZ4356O0NAXP5yGN0+X1wkGoeecMEwczWLAUAZEfVYNrM/lKS3P1RM3vLjN2xS+O5+68k/cLMXlhqWqmQ8nGTpEtKbZdIujHuc+SWWRh5njnyWK0drYsJgpKyW7MUAJAds41MX1z2/Udm7HvDAp/3zyV9xcx+LGmFpE9J+rSk15nZqEI9608v8DmA/Ko2QTAKqHOygE9Wa5YCALKjas60mf3Q3U+f+X2l20khZxqoIkMTBBuJnGkAQL3EyZn2Kt9Xug0gTUj5kZTdmqUAgOyYrc70aWb2qCSTdETpe5VuL254zwCgDrJYsxQAkB1Vg2l3P6SZHQGARunokNauTboXAIBWVEtpPAAAAAAVEEwDAAAAMRFMAwAAADERTAMAAAAxzVbNAwAAIHWKxVChZ3Q0LM7U1xdKXgJJIJgGAACZsXOntGqVNDkpTUyE2vGbNoXa8d3dSfcOeUSaBwAAyIRiMQTSxWIIpKWwjdrHx5PtH/KJYBoAAGTCwEAYka5kcjLsB5qNYBqIy10aGgrbWtoBAAcpFqVt26TNm8O2WKx+39HRqRHpmSYmwiqnQLMRTANxDQ9La9ZIGzdOBc7u4faaNWE/AKCqnTulpUulDRukyy8P26VLQ3slnZ0hR7qS9nZp+fLG9RWohmAaiKu3V1q/XurvnwqoN24Mt9evD/sBABXFyX/u65PaqkQubW1hP9BsBNNAXGbS1q1TAXVb21QgvXVr2A8AqChO/nOhEKp2FApTI9Tt7VPtHR2N6y9QDaXxgIWIAur+/qk2AmkAmFPc/OfubmlsLATbu3eH1I6+PgJpJIdgGliIKLWj3MaNBNQAMIco/7lSQD1X/nNHh7R2beP6BswHaR5AXDNzpCcnD86hBgBURP4zWgXBNBDX8PDBOdLlOdRU8wCAqsh/Rqswz/DoWVdXl4+MjCTdDeSVewiYe3unp3RUawcAHGR8nPxnZIOZ7XL3roPaCaYBAACA2VULpknzAAAAAGIimAYAAABiIpgGAAAAYiKYBgAAAGIimAYAAABiIpgGAAAAYiKYBgAAAGIimAYAAABiIpgGAAAAYiKYBgAAAGIimAYAAABiWpR0BwAA+VUsSgMD0uio1Nkp9fVJhULSvQKA2hFMAwASsXOntGqVNDkpTUxI7e3Spk3Sjh1Sd3fSvQOA2pDmAQBoumIxBNLFYgikpbCN2sfHk+0fANSKYBpA8tyloaGwraUdmTcwEEakK5mcDPsBIAsIpgEkb3hYWrNG2rhxKnB2D7fXrAn78yBHJxWjo1Mj0jNNTEi7dze3PwAQF8E00AqyHoT19krr10v9/VMB9caN4fb69WF/HuTopKKzM+RIV9LeLi1f3tz+AEBcBNNAK8h6EGYmbd06FVC3tU0F0lu3hv15kKOTir6+8GeupK0t7AeALDBP+4jVLLq6unxkZCTpbgDJmxl0bd168O0sBKTu0yOsycls9Lueyv+WkSz9DeehUjWPtjaqeQBIJzPb5e5dB7UTTAMtIutBWNb7X085OqkYHw+TDXfvDqkdfX1SR0fSvQKAg1ULpknzAFpFlCpRLiuB6MyR9cnJg9Md8iI6FuVa+Bh0dEhr10qXXRa2BNIAsoZgGmgVWQ7ChocPTkkpz6FOe853vXBSAQCZQzANtIKsB2G9vdLg4PSR9CigHhxsqYl3s+KkAgAyh5xpoBUMDYWqHeVBWHmAPTgorV6ddC8xF/cQMPf2Tk/PqdYOAGgaJiACrYwgDACAhmICItDKzMLI88yAuVo70iHri+0AAAimASAxWV9sBwCgRUk8qZntkVSU9JSkJ929y8yOkjQgaZmkPZLe6u4PJ9E/AGiK8hUPpYMX28nLxEsAyLAkR6bPdvcVZbknWyTd7O6dkm4u3QZQT6QVpAvLqANA5qUpzeNCSdeWvr9WUnqHZAhIkFWkFaRPlhfbAQAkFky7pG+a2S4zW1dqe467PyhJpe2zE+rb3AhIkFXlaQXR65e0gmRlebEdAEAyOdOSXu3uY2b2bEnfMrP7an1gKfheJ0nHHXdco/o3O/IckVXlo6D9/VOvYdIKkjHzZKb8vUTibwIAGZB4nWkz+7ikcUnvlfQad3/QzI6RdIu7v3C2xyZaZ7r8QzBCQIKscA/5uZHJSV63SWCxHQDIjNTUmTazdjMrRN9LOlfS3ZJuknRJ6W6XSLqx2X2bF/IckVXu0oYN09uitALy/puLZdQBIPOSSPN4jqQhCx8ciyT9s7t/3cz+U9JXzWytpAckvSWBvtWuWp4jATXSLHrdXnFFuP3BD4Ztf/9UAH3FFYyINku0qE6t7QCA1Gl6MO3uP5d0WoX2/ZJWNrs/sZDniKwaHg6v0yiIvuKK8P0HPzgVYJP3DwBAzZKagJhtUUBSnudYPqmrp4dRJaRTlFYQBctm0/P+P/hBTgYBNFSxKA0MSKOjUmen1NcnFQpJ9wqIL/EJiAuR2ARE9xBQ9/ZODzqqtQNpxUREAE20c6e0alV4q5mYkNrbw1vQjh1Sd3fSvQNml5oJiC0hymecGXRUawfSiPrGAJqoWAyBdLEYAmkpbKP28fFk+wfERTAN5NHMvP/JyYMXcwGAOhoYCG81lUxOhv1AFpEzDeQRef9AarVqTvHo6NSI9EwTE9Lu3c3tD1AvBNNAHpVPRJxZ37inh2oeQEIq5RRv2tQaOcWdneH3qRRQt7dLy5c3v09APTABEQCAFCgWpaVLw3amQkEaG5M6Oprfr3pp9d8PrY8JiAAApFir5xQXCmGEvVAII9FS2EbtBNLIKoLpvKi2TDTLR6MaXjNAU+Uhp7i7O4xA9/dLW7aE7dhY9lNYkG8E03kxPCytWTO9UkNU0WHNmrAfKMdrBnFwEhZblFNcSSvlFHd0SGvXSpddFraMSCPrCKbzorf34NJn5aXRmHCGmXjNII6cn4QVi9K2bdLmzWFbKT+4mr6+6WsolWtrC/sBpA8TEPOkPBiKlJdGA2biNYP5mnnStXXrwbdb9LVTj9X9WCEQSK9qExAJpvOG5aMxX7xmMF85PAmrZ6WK8fEw2XD37pDa0ddHKgSQBlTzAMtHY/54zSCO8kWAIi0cSEv1rcRBTjGQLQTTecHy0ZjLzAli5a+Z886TnnqK1wxqk8OTsDxU4gBQGcF0XlRbPjoKjlp8UhBqMHPiWPSaWbFC+trXpBtv5DWDueX0xD0vlTjQWhYyYRZTyJnOiyg4Kl8+erZ25M/MIOhv/kZ605tCIF1+EsZrBrMZGgonZTNfM9Fra3BQWr066V7WHav7IWuY7Dp/TEAEMLccThxDneX4xJ3gBFI4oRoYCKk/nZ1hAmmhkHSvpuPkLx6CaQC1oXoHEBuVOPItKydU27ZJGzZUzvNvbw/jKWvXNr9faVctmF6URGcApFS1iWOMTAM1iSpxIH+KxRBIl4/2RsHqqlXpGu1lwmx9MQERQJDTiWNA3jEJrT7qWR6x0ZgwW1+MTAMIqlV8kUJ7T09LThwD8qxSWsKmTelLS8iCLI329vWFv3MlLF0/f4xMAwh6e0OlhfKUjiigHhwM+wG0jPK0hCgInJiYah8fT7Z/WZOl0d5CIZwwFQpTfW5vn2pPSzpKVjABEQCAHGISWn1lsUIGE2bnhwmIrSzHpagAAPFkKS0hC6JR3WrVPNIYpDJhtj5I84hr5tLLc7U30syV66J+bNwY2lmpDgAwQ5bSErKiuzuMQPf3S1u2hO3YGPnnrY5gOq40BbC9vQdXXSivykCuKwBghr6+6SXlyzEJLb5otPeyy8I2jSPSqC/SPOIqD2ClMEkrqQB2ZtWFqE+sXAcAqCLNaQlZWEUQiDABcSHStvQyK9cBAOYpbZPQsrKK4HxxgpB9LCfeKPUOYONOJkxbYA8AwDxlsSJGLVr1BCFvqgXT5EwvRLWllxdyghInF5uV6wAALSBLqwjWinrerY9gOq5GBbBxJhNWW7ku+jlU8wDQCGmqaoSW0Irl+lrxBAHTEUzH1agAdubPaWs7+HlmYuU6AElIU1UjJKPOJ1StWK6vFU8QMB3BdFyNDGDLq3NEZst9NpNWrz54f7V2AI2Rt5FaynKizidUrViurxVPEDAdwXRcjQxgG5GLDaDx8jZSG+dKGlpLnU+oonJ9hcJUANrePtWexcmHrXiCgOmo5pE2M9+IZtav5gMKSK+8/v9SljPfGlBNKm3l+haKah6tgdJ4WTE0FEawyt+Iyt+oBgfDyDeAdMpbmcq8/b6ojBOqObXaCUIeEUxnRdw60wDSIy+BRV5H4jEdJ1TICepMZwWTCYFsy9OcB8pygnUOAIJpAKibvAUWlOXMXwWXmTihAkjzAIC6Yc5D/uT9b05qInKENA8AaLRWHKnN+8jrXPJea5vURIBgGgDqphUDi7zVzp4vam0DuUcwDQCoLu8jr7WY76q1AFoKwTQAoDpGXueWpwou1ZAOhBwjmAYAzI6R1+ryVsGlGtKBkGME0wCA2THyWh2l4QLSgZBji5LuAAAgxWZb5VBihDqq4FJeAi4KqHt68hNEll+96O+fen2QDoQcSKzOtJkdImlE0i/d/XwzO0HSdZKOknSHpHe6+x9n+xnUmQaABst7HWXMj3vIq49MThJIo2Wksc70ekn3lt3+jKSt7t4p6WFJaxPpFQBgSivWzkZjkA6EnEokmDazYyWdJ2lb6bZJOkfSDaW7XCuJd2gASFor1s5G/TEREzmWVM705yX9paRC6fbRkh5x9ydLt/dKWppExwAAwDxVm4gphfaeHtKB0LKaPjJtZudL2ufuu8qbK9y14mmsma0zsxEzG3nooYca0kcAGUFtWyAdSAdCjiWR5vFqSW8ysz0KEw7PURipPtLMopHyYyWNVXqwu1/p7l3u3rVkyZJm9BdAWlHbFkgH0oGQY00Ppt39I+5+rLsvk3SxpO+4+zskfVfSRaW7XSLpxmb3DUDGUNsWQF5xZS410rRoy2ZJm8xst0IO9VUJ9wdA2rHUNYC84spcaiRWZ7oeqDMNQBK1bQHkz2wLKjGg0BBprDMNAAtHbVsAecSVudQgmAaQXdS2BZBn5SUII/UKpMnJrhnBNIDsqlbbNgqoyRkE0MoaeWWOnOyaEUwDyC5q2wLIq0ZfmaNaUs2YgAgAAJA1Q0NhhLj8ylx5wDs4uPBVJ8t/XiTHOdnVJiASTAMAAGSNe0i16O2dHthWa1/I81AtSRLVPAAAAFpHM1adpFpSTQimAQAAMB3Vkmq2KOkOAAAAIGWqVUuSQntPz8JzslsEwTQAAACmi6olledeRwF1Tw/VPMoQTAMAAGC6KPe61vYcI2caAAAAiIlgGgAAAIiJYBoAAACIiWAaAAAAiIlgGgAAAIiJYBpAurlLQ0MHLxBQrR0AgCYimAaQbsPD0po101fcilbmWrMm7AcAICHUmQaQbr29U0vYSmHBgPIlblk4AACQIIJpAOk2cwnbKKguX+IWAICEmGc437Crq8tHRkaS7gaAZnCX2soy0yYnCaQBAE1jZrvcvWtmOznTANIvypEuV55DDQBAQgimAaRbFEhHOdKTk1M51ATUAICEkTMNIN2Gh6cC6ShHujyHuqdHWr062T4CAHKLYBpAuvX2SoODYRvlSEcBdU8P1TwAAIkimAaQbmaVR56rtQMA0ETkTAMAAAAxEUwDAAAAMRFMAwAAADERTAMAAAAxEUwDAAAAMRFMAwAAADERTAMAAAAxEUwDAAAAMRFMAwAAADERTAMAAAAxEUwDAAAAMRFMAwAAADERTAMAAAAxEUwDAAAAMRFMAwAAADERTAMAAAAxEUwDAAAAMRFMAwAAADERTAMAsFDu0tBQ2NbSDqBlEEwDALBQw8PSmjXSxo1TgbN7uL1mTdgPoCUtSroDAABkXm+vtH691N8fbm/dGgLp/v7Q3tubbP8ANAzBNAAAC2UWAmgpBNBRUL1+fWg3S65vABqq6WkeZrbYzP7DzH5kZj8xs78utZ9gZreb2aiZDZjZYc3uGwAAsZUH1BECaaDlJZEz/QdJ57j7aZJWSHqDmb1S0mckbXX3TkkPS1qbQN8AAIgnypEuV55DDaAlNT2Y9mC8dPPQ0pdLOkfSDaX2ayWRYAYAyIYokI5ypCcnp3KoCaiBlpZIzrSZHSJpl6Tlkv5W0s8kPeLuT5buslfS0iT6BgDAvA0PTwXSUWpHeQ51T4+0enWyfQTQEIkE0+7+lKQVZnakpCFJL650t0qPNbN1ktZJ0nHHHdewPgIAULPeXmlwMGyjHOkooO7poZoH0MISrTPt7o9IukXSKyUdaWZRcH+spLEqj7nS3bvcvWvJkiXN6SgAALMxCyPPMycbVmsH0DKSqOaxpDQiLTM7QtJrJd0r6buSLird7RJJNza7bwAAAMB8JJHmcYyka0t5022Svuru/2pm90i6zsw+KemHkq5KoG8AAABAzZoeTLv7jyWdXqH955Je3uz+AAAAAHElmjMNAAAAZBnBNAAAABATwTQAAAAQE8E0AAAAEBPBNAAAABATwTQAAAAQE8E0AAAAEBPBNAAAABATwTQAAAAQE8E0AAAAEBPBNAAAABATwTQAAAAQE8E0AAAAEBPBNAAAABCTuXvSfYjNzB6SdH/C3XiWpN8k3IdWxbFtHI5tY3BcG4dj2zgc28bh2DZOEsf2eHdfMrMx08F0GpjZiLt3Jd2PVsSxbRyObWNwXBuHY9s4HNvG4dg2TpqOLWkeAAAAQEwE0wAAAEBMBNMLd2XSHWhhHNvG4dg2Bse1cTi2jcOxbRyObeOk5tiSMw0AAADExMg0AAAAEBPBdI3MbLGZ/YeZ/cjMfmJmf11qP8HMbjezUTMbMLPDku5rVpnZIWb2QzP719Jtjm0dmNkeM7vLzO40s5FS21Fm9q3Ssf2WmT0z6X5mkZkdaWY3mNl9ZnavmZ3JsV04M3th6fUafT1qZhs4tgtnZhtLn2F3m9n20mcb77V1YGbrS8f1J2a2odTGazYGM/sHM9tnZneXtVU8lhZcYWa7zezHZnZGs/tLMF27P0g6x91Pk7RC0hvM7JWSPiNpq7t3SnpY0toE+5h16yXdW3abY1s/Z7v7irIyQlsk3Vw6tjeXbmP++iV93d1fJOk0hdcvx3aB3P2npdfrCkkvlfSYpCFxbBfEzJZK+qCkLnc/WdIhki4WvJvjPwAABxRJREFU77ULZmYnS3qvpJcrvBecb2ad4jUb1zWS3jCjrdqxfKOkztLXOklfbFIfDyCYrpEH46Wbh5a+XNI5km4otV8rqTeB7mWemR0r6TxJ20q3TRzbRrpQ4ZhKHNtYzOzpks6SdJUkufsf3f0RcWzrbaWkn7n7/eLY1sMiSUeY2SJJT5P0oHivrYcXS/p3d3/M3Z+U9H8krRav2Vjc/VZJv53RXO1YXijpy6U47d8lHWlmxzSnpwHB9DyU0hDulLRP0rck/UzSI6V/HEnaK2lpUv3LuM9L+ktJk6XbR4tjWy8u6ZtmtsvM1pXanuPuD0pSafvsxHqXXc+X9JCkq0vpSdvMrF0c23q7WNL20vcc2wVw919K+pykBxSC6N9J2iXea+vhbklnmdnRZvY0SaskPU+8Zuup2rFcKukXZfdr+muYYHoe3P2p0mXHYxUu5by40t2a26vsM7PzJe1z913lzRXuyrGN59XufobCpbBLzeyspDvUIhZJOkPSF939dEkT4hJuXZVyd98k6fqk+9IKSjmmF0o6QdJzJbUrvC/MxHvtPLn7vQrpMt+S9HVJP5L05KwPQr0kHi8QTMdQupR7i6RXKlxOWFTadayksaT6lWGvlvQmM9sj6TqFS46fF8e2Ltx9rLTdp5B3+nJJv44ug5W2+5LrYWbtlbTX3W8v3b5BIbjm2NbPGyXd4e6/Lt3m2C7MayX9t7s/5O5PSBqU9CrxXlsX7n6Vu5/h7mcppCiMitdsPVU7lnsVrgJEmv4aJpiukZktMbMjS98fofCmdK+k70q6qHS3SyTdmEwPs8vdP+Lux7r7MoVLut9x93eIY7tgZtZuZoXoe0nnKlyOvEnhmEoc21jc/VeSfmFmLyw1rZR0jzi29fQ2TaV4SBzbhXpA0ivN7GmleSnRa5b32jows2eXtsdJWqPw2uU1Wz/VjuVNkt5VqurxSkm/i9JBmoVFW2pkZqcqJLwfonAS8lV3/4SZPV9hNPUoST+U9H+7+x+S62m2mdlrJH3Y3c/n2C5c6RgOlW4ukvTP7v6/zOxoSV+VdJzCB+xb3H3mZA/MwcxWKEyaPUzSzyW9R6X3B3FsF6SUd/oLSc9399+V2njdLpCFsq59CikIP5T0pwr5pbzXLpCZfU9hvs8Tkja5+828ZuMxs+2SXiPpWZJ+LeljkoZV4ViWTgy/oFD94zFJ73H3kab2l2AaAAAAiIc0DwAAACAmgmkAAAAgJoJpAAAAICaCaQAAACAmgmkAAAAgJoJpAGgiM3vKzO4s+2raqolm9g9mts/M7p7lPi80s1tKfbvXzK5sVv8AIIsojQcATWRm4+7ekdBznyVpXNKX3f3kKvf5hqT/7e43lm6f4u53LfB5D3H3pxbyMwAgrRiZBoCEmdkzzOyn0WqKZrbdzN5b+v6LZjZiZj8pLbgRPWaPmX3KzH5Q2n+GmX3DzH5mZn9W6Xnc/VaFZY5nc4zC8rzRY+4qPd8hZvY5M7vLzH5sZn9eal9pZj8stf+DmR1e1r+PmtlOSW8xsxPN7OtmtsvMvmdmL4p/xAAgPRYl3QEAyJkjzOzOstuXufuAmX1A0jVm1i/pme7+96X9/09pla9DJN1sZqe6+49L+37h7mea2VZJ10h6taTFkn4i6Usx+7dV0nfM7PuSvinpand/RNI6SSdIOt3dnzSzo8xscel5V7r7f5nZlyW9X9LnSz/rcXfvliQzu1nSn7n7qJm9QtL/lnROzD4CQGoQTANAc/3e3VfMbHT3b5nZWyT9raTTyna91czWKbxfHyPpJElRMH1TaXuXpA53L0oqmtnjZnZkKQieF3e/upTq8QZJF0p6n5mdJum1kv7/9u6mRacwjuP49zeNpDSLsWBJiaaUlY0XIF4CstBYKCth5wV4CRa6YyV5BZZIlhoPkc3UNGTlqYSSv8V1TXQaC0d03/X9bM7TdZ1zVqffufpf51ytqm+93du+f7WqXvbuN4Bz/AzTtwCSbAcOA7fbn38B2Pqn9yZJ08gwLUlTIMkcsAR8BhaB9SR7gIvAoap6l+Q6beR5w9e+/P7L+sb26Od7Vb0GJsCkT1Y8AAQYTrLJsO/Ap76cA95v9hIhSbPOmmlJmg7ngefAcVqI3QIs0ALphyQ7gWP/+iaSHO3XJskuYAfwilbycTbJfD+2CLwAdifZ27ufAu4Oz1lVH4HVPvJOmoPDdpI0iwzTkvR/bRt8Gu9Kkn3AGeBCVd0H7gGXq2oFeESrgZ4AD/7mwkluAg+B/UnWkyxv0uwI8DTJCnAHuFRVb4BrwBrwuB87UVVfgNO08o0ntBHx39VqnwSWe99ntBISSZp5fhpPkiRJGsmRaUmSJGkkw7QkSZI0kmFakiRJGskwLUmSJI1kmJYkSZJGMkxLkiRJIxmmJUmSpJEM05IkSdJIPwApdkXio+QQygAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 864x576 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "positive = data[data['Admitted'].isin([1])]\n",
    "negative = data[data['Admitted'].isin([0])]\n",
    "\n",
    "fig, ax = plt.subplots(figsize=(12, 8))\n",
    "ax.scatter(positive['Exam 1'],\n",
    "           positive['Exam 2'],\n",
    "           s=50,\n",
    "           c='b',\n",
    "           marker='o',\n",
    "           label='Admitted')\n",
    "ax.scatter(negative['Exam 1'],\n",
    "           negative['Exam 2'],\n",
    "           s=50,\n",
    "           c='r',\n",
    "           marker='x',\n",
    "           label='Not Admitted')\n",
    "ax.legend()\n",
    "ax.set_xlabel('Exam 1 Score')\n",
    "ax.set_ylabel('Exam 2 Score')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "看起来在两类间，有一个清晰的决策边界。现在我们需要实现逻辑回归，那样就可以训练一个模型来预测结果。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Sigmoid 函数\n",
    "$g$ 代表一个常用的逻辑函数（logistic function）为$S$形函数（Sigmoid function），公式为： \\\\[g\\left( z \\right)=\\frac{1}{1+{{e}^{-z}}}\\\\] \n",
    "合起来，我们得到逻辑回归模型的假设函数： \n",
    "\t\\\\[{{h}}\\left( x \\right)=\\frac{1}{1+{{e}^{-{{w }^{T}}x}}}\\\\] "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "def sigmoid(z):\n",
    "    return 1 / (1 + np.exp(-z))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "让我们做一个快速的检查，来确保它可以工作。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAsIAAAHSCAYAAADmLK3fAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAgAElEQVR4nO3deZxcVZ3//9cnKyA7hDUsAYIIqIARXBhZhRAx7CRB+TGKwyCiP51RNpFhURCYwRVZZEANVCcQtggBEUT2KAHZEQjRQIBhC5AgZD/fP27HdDrdSSVd3aeq7uv5eNSja7nd/c593K5+5/S550ZKCUmSJKlseuUOIEmSJOVgEZYkSVIpWYQlSZJUShZhSZIklZJFWJIkSaVkEZYkSVIp9cn1jdddd920+eab5/r2kiRJKomHHnrojZTSgPbPZyvCm2++OZMmTcr17SVJklQSETG1o+edGiFJkqRSsghLkiSplCzCkiRJKiWLsCRJkkrJIixJkqRSsghLkiSplCzCkiRJKiWLsCRJkkrJIixJkqRSsghLkiSplCzCkiRJKiWLsCRJkkrJIixJkqRSsghLkiSplJZZhCPi8oh4LSKe6OT1iIifRsTkiHgsInaqfUxJkiSptqoZEf4VMHQpr+8HDG69HQNc1PVYkiRJUvfqs6wNUkp3R8TmS9nkAOA3KaUETIyINSNiw5TSKzXKKEmSpO6U0qKP3Xl/rbW699+xnJZZhKuwMfBim8fTWp+zCEuSpLwWLIB334WZM+G992Du3I5v8+Z1/trSbsvzefPnF3na3lJa8rmlPb8in9P2+Y7KaU/p1avYB3WkFkU4Oniuwz0bEcdQTJ9g0003rcG3liRJTWf+/KK4zpwJM2as2P2FH999t/ty9u27fLfevaFfv6IQ9uoFEYvut70t7/PVfk7Eohv0/P1e9bdGQy2K8DRgkzaPBwIvd7RhSulS4FKAIUOG9PB/QyRJUo9asABeegkmT4bnn4dXXqmuzL73XnVfv18/WG01WH314uNqq8GAAbDFFks+v9pqsMoqxecsLKZ9+ix/mW1batXwalGExwPHR8QYYBfgHecHS5JUEnPnwgsvFGV3YeFdeH/KFJg9e/HtV155yZK60UbwwQ8u+Xz7Itv+cf/+ef7NahrLLMIR0QLsDqwbEdOA/wL6AqSULgYmAMOAycB7wJe6K6wkScpg1qyi1LYtuQvv//3vi8/7XGUV2Gor2GYb2H//4v6WWxYfN9qoGE2V6kQ1q0aMWsbrCfhazRJJkqSeN3PmonLbvvBOm7b4iVVrrlkU249/HEaOLO4vLLwbbLBoXqhU52oxNUKSJDWC6dMXFdz2hfe11xbfdr31inK7xx6LRnQX3tZeO09+qcYswpIkNbO//Q2uvLK4Pfvs4q9tsklRcocPX3wKw5ZbFnNwpSZnEZYkqdm89RZccw2MHg333ls8t/vucMwxsPXWRdEdNKg4cU0qMYuwJEnNYM4cuOWWovz+9rfF4w99CM4+G77wBXD9fmkJFmFJkhpVSvCnPxXld+xYePPNYm7vV78KRx4JO+3kiWvSUliEJUlqNFOmLJr3+9xzsNJKcOCBRfndZ5/iQhGSlsmfFEmSGsFbb8HVVxejv/fdV4z07r47nHwyHHJIcbEJScvFIixJUr2aMwcmTCjK7003FY+33RbOOaeY97vJJrkTSg3NIixJUj1JCSZOXDTvd/r0Yt7vcccVUx923NF5v1KNWIQlSaoHzz+/aN7v5MnF0mYL5/1+9rPO+5W6gT9VkiTlsnDe729+A/ffX4z07rEHfPe7cPDBzvuVuplFWJKknjR7djHv98orF8373W47+OEP4YgjnPcr9SCLsCRJPeGpp+DnP18073f99eFrXyumPuywg/N+pQwswpIkdbebboIRI4oT4Q46qCi/e+/tvF8pM38CJUnqTr/8JRx7bLHaw003wQYb5E4kqVWv3AEkSWpKKcFpp8Exx8C++8If/2gJluqMI8KSJNXa3LlFAf7Vr+DLX4aLL4a+fXOnktSOI8KSJNXSzJnw+c8XJfj00+GyyyzBUp1yRFiSpFr5v/+Dz30OHn20KMBHH507kaSlsAhLklQLzzwDQ4fCa6/B+PEwbFjuRJKWwSIsSVJX3XcfDB9eLId2110wZEjuRJKq4BxhSZK64vrrizWB11kHHnjAEiw1EIuwJEkr6uc/h0MOKa4Md//9sMUWuRNJWg4WYUmSlteCBXDiifD1rxdTIu64A9ZdN3cqScvJOcKSJC2P2bOLtYErFfjqV+FnP4PevXOnkrQCLMKSJFXrnXfgoIPgzjvhnHOKUeGI3KkkrSCLsCRJ1Zg2rVgS7emn4Te/gSOPzJ1IUhdZhCVJWpYnnoD99itGhG+5pVglQlLD82Q5SZKW5o9/hF13hfnz4e67LcFSE7EIS5LUmTFjYN99YaONYOLEYpk0SU3DIixJUnspwf/8D4waBbvsUlw5btNNc6eSVGMWYUmS2po/H771Lfj2t+Gww+C222CttXKnktQNPFlOkqSFZs2CL34Rrr0WvvnNYlS4l2NGUrOyCEuSBDB9OhxwANx7L1xwQTEqLKmpWYQlSZo6tVge7fnnYexYOPzw3Ikk9QCLsCSp3B55pCjBs2YV84F32y13Ikk9xIlPkqTy+v3v4V/+Bfr0KaZEWIKlUrEIS5LK6Te/KS6ZPGhQsUbwdtvlTiSph1mEJUnlkhKcfTYcdRR85jNwzz2w8ca5U0nKwDnCkqTymDcPvv51uPhiOOIIuOIK6NcvdypJmTgiLEkqh/feg4MPLkrwiSfC6NGWYKnkHBGWJDW/lIqrxN1yC/z85/C1r+VOJKkOWIQlSc3vtttgwgT47/+2BEv6J6dGSJKa24IFxVSIQYOK+cGS1MoRYUlSc6tU4NFHi4/OCZbUhiPCkqTmNWsWnHoq7LQTjBiRO42kOuOIsCSpeV10EUydCpddBr0c+5G0ON8VJEnN6e234fvfh332gb33zp1GUh2yCEuSmtN558H06XDuubmTSKpTFmFJUvN56SX48Y/hC1+AHXbInUZSnbIIS5Kaz+mnw/z5xdQISeqERViS1FyeegouvxyOOw423zx3Gkl1zCIsSWoup5wCq64K3/1u7iSS6pxFWJLUPO67D268sbiS3Lrr5k4jqc5ZhCVJzSElOOEE2HBD+OY3c6eR1AC8oIYkqTnceCPcfz9ceimsskruNJIagCPCkqTGN28enHwybLMNfOlLudNIahCOCEuSGt8VV8Bf/wrXXw99/NUmqTqOCEuSGts//gH/9V/wqU/BAQfkTiOpgfjfZklSY/vJT+CVV+CaayAidxpJDcQRYUlS43rjDTj33GIk+NOfzp1GUoOxCEuSGtcPfgDvvgtnn507iaQGZBGWJDWmv/0NLrwQvvxl2Hbb3GkkNSCLsCSpMX3ve8UKEaefnjuJpAZlEZYkNZ6//AWuuqq4gtzGG+dOI6lBWYQlSY3npJNg7bWLSypL0gpy+TRJUmO5/Xa47Ta44AJYc83caSQ1MEeEJUmNY8ECOPFE2GwzOO643GkkNThHhCVJjWPsWHj4YRg9Gvr3z51GUoNzRFiS1BjmzIHvfhc++lE44ojcaSQ1gaqKcEQMjYhnImJyRJzUweubRsSdEfGXiHgsIobVPqokqdQuvrhYO/jcc6GX4ziSum6Z7yQR0Ru4ENgP2BYYFRHtVy4/Fbg6pbQjMBL4Ra2DSpJKbMYMOOss2HNP2Gef3GkkNYlq/ku9MzA5pTQlpTQHGAMc0G6bBKzeen8N4OXaRZQkld7558MbbxSjwRG500hqEtWcLLcx8GKbx9OAXdptczpwW0R8HfgAsHdN0kmS9MorxVJpI0bAkCG500hqItWMCHf0X+/U7vEo4FcppYHAMGB0RCzxtSPimIiYFBGTXn/99eVPK0kqnzPOgLlz4Qc/yJ1EUpOppghPAzZp83ggS059OBq4GiCl9ACwErBu+y+UUro0pTQkpTRkwIABK5ZYklQezzwDl10Gxx4LW26ZO42kJlNNEX4QGBwRgyKiH8XJcOPbbfMCsBdARHyIogg75CtJ6ppTToGVV4ZTT82dRFITWmYRTinNA44Hfgc8TbE6xJMRcWZEDG/d7D+Bf4uIR4EW4F9TSu2nT0iSVL2JE+G66+CEE2C99XKnkdSEIldfHTJkSJo0aVKW7y1JqnMpwW67wbPPwuTJsOqquRNJamAR8VBKaYmzbb3EsiSp/tx8M9xzD1x0kSVYUrfx0jySpPoyfz6cdBJsvTUcfXTuNJKamCPCkqT68utfw5NPwrhx0Ldv7jSSmpgjwpKk+vH++3DaabDLLnDwwbnTSGpyjghLkurHT38KL70EV13lpZQldTtHhCVJ9WH6dDjnHPjc54oVIySpm1mEJUn14eyzYcYM+OEPcyeRVBIWYUlSflOnws9+BkcdBdtvnzuNpJKwCEuS8jvttGJO8Jln5k4iqUQswpKkvB57DEaPhm98AzbZJHcaSSViEZYk5XXSSbDGGnDyybmTSCoZl0+TJOVz551wyy1w3nmw1lq500gqGUeEJUl5pAQnnlhMh/j613OnkVRCjghLkvIYNw4efBCuuAJWWil3Gkkl5IiwJKnnzZ0Lp5xSLJV25JG500gqKUeEJUk975e/hMmT4aaboHfv3GkklZQjwpKknjVzJpxxRnEZ5WHDcqeRVGKOCEuSetYFF8Brr8H48cVFNCQpE0eEJUk959VX4fzz4dBDYZddcqeRVHIWYUlSzznrLJg1C37wg9xJJMkiLEnqIc89B5dcAsccA1tvnTuNJFmEJUk95LTToH//4qMk1QGLsCSp+02fDtdeW4wGb7BB7jSSBFiEJUk94dpri4tofPGLuZNI0j9ZhCVJ3a9SgQ9+EHbcMXcSSfoni7AkqXu99BLcdRcccYTrBkuqKxZhSVL3GjsWUoJRo3InkaTFWIQlSd2rUoGPfxwGD86dRJIWYxGWJHWfZ56Bhx5yNFhSXbIIS5K6T0tLMS94xIjcSSRpCRZhSVL3SKmYFrHHHrDRRrnTSNISLMKSpO7x8MPFZZWPOCJ3EknqkEVYktQ9KhXo2xcOPjh3EknqkEVYklR78+fDmDEwbBistVbuNJLUIYuwJKn27r4bXn7ZaRGS6ppFWJJUey0tsOqqsP/+uZNIUqcswpKk2po9G8aNgwMPhFVWyZ1GkjplEZYk1dbvfgdvveW0CEl1zyIsSaqtSgXWXRf23jt3EklaKouwJKl23n0Xxo+Hww8vlk6TpDpmEZYk1c6NN8L77zstQlJDsAhLkmqnUoFNN4VPfjJ3EklaJouwJKk2Xn+9OFFu1Cjo5a8XSfXPdypJUm2MG1dcUc5pEZIahEVYklQbLS2w3Xbw4Q/nTiJJVbEIS5K67oUX4J57imkREbnTSFJVLMKSpK4bM6b4OGpU3hyStBwswpKkrqtU4BOfgC22yJ1EkqpmEZYkdc1TT8Gjj3qSnKSGYxGWJHVNS0uxXNrhh+dOIknLxSIsSVpxKRXTIvbaC9ZfP3caSVouFmFJ0or7859hyhSnRUhqSBZhSdKKq1Sgf3846KDcSSRpuVmEJUkrZv58GDsW9t8f1lgjdxpJWm4WYUnSirnzTnj1VdcOltSwLMKSpBVTqcDqq8OwYbmTSNIKsQhLkpbfrFlw7bVw8MGw8sq500jSCrEIS5KW3y23wIwZrhYhqaFZhCVJy69SgfXWgz32yJ1EklaYRViStHxmzIDf/hZGjIA+fXKnkaQVZhGWJC2f66+H2bOdFiGp4VmEJUnLp1KBQYNgl11yJ5GkLrEIS5Kq9+qrcMcdxWhwRO40ktQlFmFJUvWuuaa4opwX0ZDUBCzCkqTqVSrwkY/AdtvlTiJJXWYRliRV529/gwce8CQ5SU3DIixJqs6YMcXHkSPz5pCkGrEIS5KqU6nApz8Nm22WO4kk1YRFWJK0bI8/Dk884bQISU3FIixJWrZKBXr3hsMOy51EkmqmqiIcEUMj4pmImBwRJ3WyzeER8VREPBkRldrGlCRlkxK0tMA++8CAAbnTSFLNLPMi8RHRG7gQ+CwwDXgwIsanlJ5qs81g4GTg0ymltyJive4KLEnqYQ88AFOnwve/nzuJJNVUNSPCOwOTU0pTUkpzgDHAAe22+TfgwpTSWwAppddqG1OSlE2lAiutBAe0f+uXpMZWTRHeGHixzeNprc+1tTWwdUTcFxETI2JorQJKkjKaOxeuvhqGD4fVVsudRpJqaplTI4COLiafOvg6g4HdgYHAPRGxfUrp7cW+UMQxwDEAm2666XKHlST1sDvugNdfd7UISU2pmhHhacAmbR4PBF7uYJsbU0pzU0p/A56hKMaLSSldmlIaklIaMsATLiSp/rW0wJprwlD/0Cep+VRThB8EBkfEoIjoB4wExrfb5gZgD4CIWJdiqsSUWgaVJPWw99+H666DQw6B/v1zp5GkmltmEU4pzQOOB34HPA1cnVJ6MiLOjIjhrZv9DngzIp4C7gS+k1J6s7tCS5J6wE03wbvvOi1CUtOKlNpP9+0ZQ4YMSZMmTcryvSVJVTjoIPjTn+DFF4uLaUhSg4qIh1JKQ9o/75XlJElLevttmDABRo60BEtqWhZhSdKSrrsO5syBUaNyJ5GkbmMRliQtqVKBrbaCIUv8JVGSmoZFWJK0uFdegT/8oThJLjpaSl6SmoNFWJK0uLFjISWnRUhqehZhSdLiWlpgp51gm21yJ5GkbmURliQtMnky/PnPjgZLKgWLsCRpkZaWYl7wyJG5k0hSt7MIS5IKKRWrRXzmMzBwYO40ktTtLMKSpMKjj8Jf/+ollSWVhkVYklSoVKBPHzjkkNxJJKlHWIQlSbBgQTE/eOhQWGed3GkkqUdYhCVJcO+9MG2a0yIklYpFWJJUjAavsgoMH547iST1GIuwJJXdnDlw9dVw4IHwgQ/kTiNJPcYiLEll9/vfw/TpXkRDUulYhCWp7CoVWHtt2Gef3EkkqUdZhCWpzP7xD7jhBjjsMOjXL3caSepRFmFJKrPf/hbee8/VIiSVkkVYksqsUikup7zrrrmTSFKPswhLUlm9+SbccguMHAm9/HUgqXx855Oksrr2Wpg3z2kRkkrLIixJZdXSAttsAzvskDuJJGVhEZakMpo2De66q1g7OCJ3GknKwiIsSWU0diyk5EU0JJWaRViSyqhSgY9/HAYPzp1EkrKxCEtS2TzzDDz8sCfJSSo9i7AklU1LSzEveMSI3EkkKSuLsCSVSUrFtIg99oANN8ydRpKysghLUpk89BA895zTIiQJi7AklUulAv36wcEH504iSdlZhCWpLObPL5ZNGzYM1lordxpJys4iLEllcffd8PLLrh0sSa0swpJUFpUKrLoq7L9/7iSSVBcswpJUBrNnw7hxcOCBsMoqudNIUl2wCEtSGdx6K7z9NnzhC7mTSFLdsAhLUhm0tMCAAbDXXrmTSFLdsAhLUrObORPGj4fDDoO+fXOnkaS6YRGWpGZ3443w/vteREOS2rEIS1Kzq1Rgs83gk5/MnUSS6opFWJKa2euvw223FWsH9/ItX5La8l1RkprZuHHFFeW8iIYkLcEiLEnNrFKB7baDD384dxJJqjsWYUlqVi+8APfeW5wkF5E7jSTVHYuwJDWrMWOKj06LkKQOWYQlqVlVKvCJT8CgQbmTSFJdsghLUjN68kl49FHXDpakpbAIS1Izamkplks7/PDcSSSpblmEJanZpFQU4b33hvXXz51GkuqWRViSms2f/wxTpjgtQpKWwSIsSc2mUoH+/eGgg3InkaS6ZhGWpGYybx6MHQv77w+rr547jSTVNYuwJDWTO++EV191WoQkVcEiLEnNpKWlGAkeNix3EkmqexZhSWoWs2bBtdfCwQfDSivlTiNJdc8iLEnNYsIEmDHDaRGSVCWLsCQ1i0qlWDd4jz1yJ5GkhmARlqRm8M47cNNNMGIE9OmTO40kNQSLsCQ1gxtugNmznRYhScvBIixJzaBSgS22gJ13zp1EkhqGRViSGt2rr8Ltt8OoURCRO40kNQyLsCQ1uquvhgULnBYhScvJIixJja6lBT76Udh229xJJKmhWIQlqZFNmQIPPFBMi5AkLReLsCQ1sjFjio8jR+bNIUkNyCIsSY0qJbjqKth1V9hss9xpJKnhWIQlqVE9/jg89ZQnyUnSCrIIS1KjammB3r3h0ENzJ5GkhmQRlqRGtGBBUYT32QcGDMidRpIakkVYkhrRAw/A1KlOi5CkLqiqCEfE0Ih4JiImR8RJS9nu0IhIETGkdhElSUuoVGDlleGAA3InkaSGtcwiHBG9gQuB/YBtgVERscSq7RGxGvAN4E+1DilJamPuXLjmGhg+HFZbLXcaSWpY1YwI7wxMTilNSSnNAcYAHQ1BnAWcB8yqYT5JUnt33AGvv+5FNCSpi6opwhsDL7Z5PK31uX+KiB2BTVJKN9UwmySpI5UKrLkmDB2aO4kkNbRqinB08Fz654sRvYAfAf+5zC8UcUxETIqISa+//nr1KSVJhffeg+uvL5ZM698/dxpJamjVFOFpwCZtHg8EXm7zeDVge+CPEfF34BPA+I5OmEspXZpSGpJSGjLA5X4kafnddBO8+66rRUhSDVRThB8EBkfEoIjoB4wExi98MaX0Tkpp3ZTS5imlzYGJwPCU0qRuSSxJZdbSAhtuCJ/5TO4kktTwllmEU0rzgOOB3wFPA1enlJ6MiDMjYnh3B5QktXrrLZgwAUaOLK4oJ0nqkj7VbJRSmgBMaPfcaZ1su3vXY0mSlnDddTBnjtMiJKlGvLKcJDWKSgUGD4aPfSx3EklqChZhSWoEr7wCd95ZjAZHR4v5SJKWl0VYkhrB2LGQkhfRkKQasghLUiOoVGCnneCDH8ydRJKahkVYkurdc8/Bgw96kpwk1ZhFWJLqXUtLMS94xIjcSSSpqViEJamepVRMi/jMZ2DgwNxpJKmpWIQlqZ498gg884zTIiSpG1iEJameVSrQty8cckjuJJLUdCzCklSvFiyAMWNg6FBYZ53caSSp6ViEJale3XsvTJvmtAhJ6iYWYUmqV5UKrLIKfP7zuZNIUlOyCEtSPZozB665Bg48ED7wgdxpJKkpWYQlqR7ddhtMn+60CEnqRhZhSapHlUpxgtw+++ROIklNyyIsSfXmH/+AG2+EQw8tlk6TJHULi7Ak1Zvx4+G995wWIUndzCIsSfWmUikup7zrrrmTSFJTswhLUj1580249VYYNQp6+RYtSd3Jd1lJqifXXgvz5hVFWJLUrSzCklRPKhXYZhvYYYfcSSSp6VmEJalevPgi3H13cZJcRO40ktT0LMKSVC/GjoWUnBYhST3EIixJ9aJSgZ13hq22yp1EkkrBIixJ9eCvf4W//MXRYEnqQRZhSaoHLS3FvOARI3InkaTSsAhLUm4pFdMi9twTNtwwdxpJKg2LsCTlNmkSTJ7sJZUlqYdZhCUpt5YW6NcPDj44dxJJKhWLsCTlNH8+jBkDw4bBmmvmTiNJpWIRlqSc7roLXnnFaRGSlIFFWJJyqlRg1VVh//1zJ5Gk0rEIS1Ius2fDuHHF3OCVV86dRpJKxyIsSbnceiu8844X0ZCkTCzCkpRLpQIDBsBee+VOIkmlZBGWpBxmzoTx4+Hww6Fv39xpJKmULMKSlMMNN8CsWa4WIUkZWYQlKYeWFthsM/jkJ3MnkaTSsghLUk97/XW47bbiJLmI3GkkqbQswpLU0665priinNMiJCkri7Ak9bRKBbbfHj784dxJJKnULMKS1JOmToX77nPtYEmqAxZhSepJY8YUHy3CkpSdRViSelKlUqwUMWhQ7iSSVHoWYUnqKU88AY895klyklQnLMKS1FNaWqB3bzjssNxJJElYhCWpZ8yZA1deCXvtBeuvnzuNJAnokzuAJJXCxRfDCy/AJZfkTiJJauWIsCR1txkz4KyzYM89Yd99c6eRJLWyCEtSdzv/fHjjDTj3XC+pLEl1xCIsSd3plVfgggtgxAgYMiR3GklSGxZhSepOZ5xRnCj3/e/nTiJJasciLEnd5Zln4LLL4NhjYautcqeRJLVjEZak7nLKKbDyyvC97+VOIknqgEVYkrrDxIlw3XXwne/AeuvlTiNJ6oBFWJJqLSU44YTiwhn/8R+500iSOuEFNSSp1m6+Ge65B37xC1h11dxpJEmdcERYkmpp/nw46SQYPBi+8pXcaSRJS+GIsCTV0q9/DU8+CddcA3375k4jSVoKR4QlqVbefx9OOw123hkOOSR3GknSMjgiLEm18tOfwksvwVVXeSllSWoAjghLUi1Mnw7nnAOf+xzstlvuNJKkKliEJakWzj4bZswoyrAkqSFYhCWpq6ZOhZ/9DI46Cj784dxpJElVsghLUleddloxJ/iMM3InkSQtB4uwJHXFY4/B6NHwjW/AppvmTiNJWg4WYUnqipNOgjXWKD5KkhqKy6dJ0oq680645RY47zxYe+3caSRJy8kRYUlaESnBiSfCwIFw/PG500iSVoAjwpK0IsaNgwcfhCuugJVXzp1GkrQCqhoRjoihEfFMREyOiCUmwkXEf0TEUxHxWETcERGb1T6qJNWJuXPhlFNg++3hyCNzp5EkraBlFuGI6A1cCOwHbAuMioht2232F2BISukjwDjgvFoHlaS68ctfwuTJ8MMfQu/eudNIklZQNSPCOwOTU0pTUkpzgDHAAW03SCndmVJ6r/XhRGBgbWNKUp2YObNYL3i33WDYsNxpJEldUM0c4Y2BF9s8ngbsspTtjwZu6UooSapbF1wAr70G48cXF9GQJDWsaopwR+/0qcMNI74IDAF26+T1Y4BjADZ14XlJjebVV+H88+HQQ2GXpY0HSJIaQTVTI6YBm7R5PBB4uf1GEbE38F1geEppdkdfKKV0aUppSEppyIABA1YkryTlc9ZZMGsW/OAHuZNIkmqgmiL8IDA4IgZFRD9gJDC+7QYRsSNwCUUJfq32MSUps+eeg0sugWOOga23zp1GklQDyyzCKaV5wPHA74CngatTSk9GxJkRMbx1s/OBVYFrIuKRiBjfyZeTpMZ06qnQvz+cdlruJJKkGqnqghoppQnAhCZP0J4AAA/gSURBVHbPndbm/t41ziVJ9ePBB+Hqq4sSvMEGudNIkmrESyxL0tKkBCecAAMGwLe/nTuNJKmGvMSyJC3NrbfCH/8IP/sZrLZa7jSSpBpyRFiSOjN/Ppx4Imy5ZXGSnCSpqTgiLEmdueoqePxxGDMG+vXLnUaSVGOOCEtSR2bNgu99Dz72MTjssNxpJEndwBFhSerIhRfCCy/AFVdAL8cMJKkZ+e4uSe299VZx9bh994U998ydRpLUTSzCktTeuefC228XHyVJTcsiLEltvfgi/OQn8MUvwkc/mjuNJKkbWYQlqa3TT4cFC+DMM3MnkSR1M4uwJC305JPwq1/B8cfD5pvnTiNJ6mYWYUla6OSTi6vHnXJK7iSSpB5gEZYkgHvugd/+Fk46CdZZJ3caSVIPsAhLUkpwwgmw8cbwjW/kTiNJ6iFeUEOSrr8eJk6Eyy6DVVbJnUaS1EMcEZZUbvPmFXODP/QhOOqo3GkkST3IEWFJ5fa//wvPPgs33gh9fEuUpDJxRFhSef3jH8W6wZ/+NHz+87nTSJJ6mMMfksrrRz+C//s/uPZaiMidRpLUwxwRllROr78O550HBx4In/pU7jSSpAwswpLK6fvfL6ZGnHNO7iSSpEwswpLKZ8oUuOgiOPpo2Gab3GkkSZlYhCWVz6mnFitEnH567iSSpIwswpLK5eGHoaUFvvUt2Gij3GkkSRlZhCWVy4knwjrrFJdUliSVmsunSSqP3/8ebr8dfvxjWGON3GkkSZk5IiypHN54o5gOsfnmcOyxudNIkuqAI8KSmt+UKbDffjB1anEp5f79cyeSJNUBi7Ck5vbQQzBsGMydC3fcUVxOWZIknBohqZndcgvsthusvDLcf78lWJK0GIuwpOZ0+eXw+c/D1lvDAw944QxJ0hIswpKaS0pwxhnFVeP22gvuugs23DB3KklSHXKOsKTmMW8efPWrcNllcNRR8MtfQt++uVNJkuqUI8KSmsO778IBBxQl+NRT4YorLMGSpKVyRFhS43v1Vdh//+LyyRdfDP/+77kTSZIagEVYUmN79tlijeBXXoEbbihOkJMkqQoWYUmNa+LEYiQ4Av74R9h559yJJEkNxDnCkhrT+PGw556w5prF8miWYEnScrIIS2o8F10EBx0E229fXChjq61yJ5IkNSCLsKTGkRKccgocd1wxL/jOO2G99XKnkiQ1KOcIS2oMc+bAV74Co0fDv/0b/OIX0Me3MEnSinNEWFL9mzEDPve5ogSfeSZccoklWJLUZf4mkVTfXn4Zhg2DJ56Ayy+HL30pdyJJUpOwCEuqX089VcwFfvNNuPlm2Hff3IkkSU3EIiypPt1zDwwfDv37w913w0475U4kSWoyzhGWVH/GjYPPfhbWX79YI9gSLEnqBhZhSfXlxz+Gww+Hj30M7rsPBg3KnUiS1KQswpLqw4IF8O1vw7e+BQceCLffDuuskzuVJKmJOUdYUn6zZ8O//iuMGQPHH1+MCvfunTuVJKnJWYQl5fX228UI8F13wbnnwne+AxG5U0mSSsAiLCmfF18slkd79lm46io44ojciSRJJWIRlpTH448XJXjmTLj1Vthzz9yJJEkl48lyknreH/4Au+5a3L/nHkuwJCkLi7CknlWpwNChsMkmxRrBH/lI7kSSpJKyCEvqGW+9BaefDl/4AnzqU3DvvUUZliQpE+cIS+o+c+bAhAkwejTcdFPxeORI+NWviksnS5KUkUVYUm2lBBMnFuV37FiYPh3WWw+++lU48sjicskujyZJqgMWYUm18fzzcOWVxW3yZFhppWJ94COPhH32gT6+3UiS6ou/mSStuOnT4eqri9Hf++8vRnp33x1OOQUOOQRWXz13QkmSOmURlrR8Zs9eNO/35puLeb/bbgvnnFOcCOcJcJKkBmERlrRsKRVLnS2c9/vWW7D++nDcccXUhx13dN6vJKnhWIQlde7554vye+WVxf2VV1407/ezn3XerySpoflbTNLipk8vRn1Hjy5GgSNgjz3g1FPh4IOd9ytJahoWYUnFvN+bb14073fuXNhuO/jhD4t5vwMH5k4oSVLNWYSlskqpWOlh9Ohi5YeF836PP76Y+rDDDs77lSQ1NYuwVDaTJy+a9ztlSjHv96CDivK7997O+5UklYa/8aRmlBK89lpReidPLk50mzwZnn4aHnmkGOndc0847bRi3u9qq+VOLElSj7MIS41qwQKYNm1RyW1beJ9/Ht59d9G2vXrBZpvBVlvBuefCEUc471eSVHoWYamezZ0LU6d2XHanTClOcluoXz/YYgvYcsvi6m5bbVXc32qrogT365ftnyFJUj2yCEu5zZpVlNq2JXfh/b//HebPX7TtKqsUxXabbWD//RcvuwMHQu/e2f4ZkiQ1mqqKcEQMBX4C9AYuSyn9sN3r/YHfAB8D3gRGpJT+XtuoUgOZPRtmzoQZM4qPC++/886SI7wvvVTM6V1ojTVg8GAYMgRGjixK7sLCu8EGruQgSVKNLLMIR0Rv4ELgs8A04MGIGJ9SeqrNZkcDb6WUtoqIkcC5wIjuCCx1i5SKkdn2xbWjMlvN/blzl/791luvKLd77LGo6C4su2uvbdmVJKkHVDMivDMwOaU0BSAixgAHAG2L8AHA6a33xwE/j4hIqe0wl0ohpeIkrgULFr/f9jZvXlEUl+e2Ip/T2ee++27HJbbtFISlWXXV4upqq61W3FZfvZibu/B+2+c7em6TTVylQZKkOlBNEd4YeLHN42nALp1tk1KaFxHvAOsAb9QiZM088QR8+9uLHrfv6W0fr+hr1W6bUvff7+z19gW1s8La0fPLei63Pn2gb9/i1vZ+29uqqxZFdKONliyty7r/gQ8UKzBIkqSGV00R7uhvtO1HeqvZhog4BjgGYNNNN63iW9fYvHnw9tvtQ3X+eEVfW9q2EYse98T9jp7r1WvRrf3jpT3flefaPt9ZQe3sVu32ffo4pUCSJFWtmiI8DdikzeOBwMudbDMtIvoAawDT23+hlNKlwKUAQ4YM6flpEzvsABMn9vi3lSRJUv2p5m+8DwKDI2JQRPQDRgLj220zHjiq9f6hwB+cHyxJkqR6tswR4dY5v8cDv6NYPu3ylNKTEXEmMCmlNB74X2B0REymGAke2Z2hJUmSpK6qah3hlNIEYEK7505rc38WcFhto0mSJEndx9PfJUmSVEoWYUmSJJWSRViSJEmlZBGWJElSKVmEJUmSVEoWYUmSJJWSRViSJEmlZBGWJElSKVmEJUmSVEoWYUmSJJWSRViSJEmlZBGWJElSKVmEJUmSVEoWYUmSJJWSRViSJEmlFCmlPN844nVgapZvDusCb2T63s3A/dc17r+ucf91jfuva9x/XeP+6xr334rbLKU0oP2T2YpwThExKaU0JHeORuX+6xr3X9e4/7rG/dc17r+ucf91jfuv9pwaIUmSpFKyCEuSJKmUylqEL80doMG5/7rG/dc17r+ucf91jfuva9x/XeP+q7FSzhGWJEmSyjoiLEmSpJJryiIcEYdFxJMRsSAihrR77eSImBwRz0TEvp18/qCI+FNEPBcRYyOiX88kr0+t++CR1tvfI+KRTrb7e0Q83rrdpJ7OWa8i4vSIeKnNPhzWyXZDW4/LyRFxUk/nrFcRcX5E/DUiHouI6yNizU628/hrY1nHU0T0b/3Zntz6frd5z6esTxGxSUTcGRFPt/4u+f872Gb3iHinzc/1aTmy1qtl/TxG4aetx99jEbFTjpz1KCI+2Oa4eiQiZkTEN9tt4/FXI31yB+gmTwAHA5e0fTIitgVGAtsBGwG3R8TWKaX57T7/XOBHKaUxEXExcDRwUffHrk8ppREL70fE/wDvLGXzPVJKrnG4pB+llP67sxcjojdwIfBZYBrwYESMTyk91VMB69jvgZNTSvMi4lzgZODETrb1+KPq4+lo4K2U0lYRMZLifW/Ekl+tlOYB/5lSejgiVgMeiojfd/DzeE9Kaf8M+RrF0n4e9wMGt952ofgdu0tPBatnKaVngB3gnz/LLwHXd7Cpx18NNOWIcErp6dYDqb0DgDEppdkppb8Bk4Gd224QEQHsCYxrferXwIHdmbdRtO6bw4GW3Fma0M7A5JTSlJTSHGAMxfFaeiml21JK81ofTgQG5szTIKo5ng6geH+D4v1ur9af8dJLKb2SUnq49f5M4Glg47ypms4BwG9SYSKwZkRsmDtUHdoLeD6llOsCZE2vKYvwUmwMvNjm8TSWfHNbB3i7zS/ejrYpq38BXk0pPdfJ6wm4LSIeiohjejBXIzi+9c9/l0fEWh28Xs2xKfgycEsnr3n8LVLN8fTPbVrf796heP9TG61TRnYE/tTBy5+MiEcj4paI2K5Hg9W/Zf08+p5XnZF0Pvjk8VcDDTs1IiJuBzbo4KXvppRu7OzTOniu/bIZ1WzTdKrcn6NY+mjwp1NKL0fEesDvI+KvKaW7a521Hi1t/1H8ye8siuPoLOB/KArdYl+ig89t+uNuoWqOv4j4LsWfrK/q5MuU9vjrgO91NRARqwLXAt9MKc1o9/LDFJdsfbd13v8NFH/mV2FZP48ef8vQen7ScIrpYO15/NVIwxbhlNLeK/Bp04BN2jweCLzcbps3KP5E06d1lKSjbZrOsvZnRPShmHf9saV8jZdbP74WEddT/Hm2FEWk2uMxIn4J3NTBS9Ucm02riuPvKGB/YK/UyZqPZT7+OlDN8bRwm2mtP99rANN7Jl79i4i+FCX4qpTSde1fb1uMU0oTIuIXEbGuc9QLVfw8lvo9r0r7AQ+nlF5t/4LHX+2UbWrEeGBk69nSgyj+9/Tnthu0/pK9Ezi09amjgM5GmMtkb+CvKaVpHb0YER9oPamEiPgAsA/FSYul127e20F0vF8eBAZHsWJJP4o/h43viXz1LiKGUpwcNzyl9F4n23j8La6a42k8xfsbFO93f+jsPxll0zpX+n+Bp1NKF3SyzQYL51RHxM4Uv0/f7LmU9avKn8fxwP/XunrEJ4B3Ukqv9HDUetfpX2E9/mqnYUeElyYiDgJ+BgwAbo6IR1JK+6aUnoyIq4GnKP7E+rWFK0ZExATgK63/iz0RGBMR3wf+QvGGWHZLzFOKiI2Ay1JKw4D1getbfy77AJWU0q09nrI+nRcRO1D82e/vwL/D4vuvdUWE44HfAb2By1NKT+YKXGd+DvSn+PMqwMSU0rEef53r7HiKiDOBSSml8RTva6MjYjLFSPDIfInrzqeBI4HHY9FykacAmwKklC6m+M/DVyNiHvA+MNL/SPxThz+PEXEs/HP/TQCGUZy0/h7wpUxZ61JErEKx6su/t3mu7f7z+KsRrywnSZKkUirb1AhJkiQJsAhLkiSppCzCkiRJKiWLsCRJkkrJIixJkqRSsghLkiSplCzCkiRJKiWLsCRJkkrp/wGXgHu7s0Y44wAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 864x576 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "nums = np.arange(-10, 10, step=1)\n",
    "\n",
    "fig, ax = plt.subplots(figsize=(12, 8))\n",
    "ax.plot(nums, sigmoid(nums), 'r')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "棒极了！现在，我们需要编写代价函数来评估结果。\n",
    "代价函数：\n",
    "$J\\left(w\\right)=-\\frac{1}{m}\\sum\\limits_{i=1}^{m}{({{y}^{(i)}}\\log \\left( {h}\\left( {{x}^{(i)}} \\right) \\right)+\\left( 1-{{y}^{(i)}} \\right)\\log \\left( 1-{h}\\left( {{x}^{(i)}} \\right) \\right))}$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "def cost(w, X, y):\n",
    "    w = np.matrix(w)\n",
    "    X = np.matrix(X)\n",
    "    y = np.matrix(y)\n",
    "    first = np.multiply(-y, np.log(sigmoid(X * w.T)))\n",
    "    second = np.multiply((1 - y), np.log(1 - sigmoid(X * w.T)))\n",
    "    return np.sum(first - second) / (len(X))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "现在，我们要做一些设置，和我们在练习1在线性回归的练习很相似。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "# add a ones column - this makes the matrix multiplication work out easier\n",
    "data.insert(0, 'Ones', 1)\n",
    "\n",
    "# set X (training data) and y (target variable)\n",
    "cols = data.shape[1]\n",
    "X = data.iloc[:, 0:cols - 1]\n",
    "y = data.iloc[:, cols - 1:cols]\n",
    "\n",
    "# convert to numpy arrays and initalize the parameter array w\n",
    "X = np.array(X.values)\n",
    "y = np.array(y.values)\n",
    "w = np.zeros(3)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "让我们来检查矩阵的维度来确保一切良好。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "((100, 3), (3,), (100, 1))"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "X.shape, w.shape, y.shape"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "让我们计算初始化参数的代价函数($w$为0)。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.6931471805599453"
      ]
     },
     "execution_count": 9,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "cost(w, X, y)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "看起来不错，接下来，我们需要一个函数来计算我们的训练数据、标签和一些参数$w$的梯度。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# gradient descent(梯度下降)\n",
    "* 这是批量梯度下降（batch gradient descent）  \n",
    "* 转化为向量化计算： $\\frac{1}{m} X^T( Sigmoid(XW) - y )$\n",
    "$$\\frac{\\partial J\\left( w  \\right)}{\\partial {{w }_{j}}}=\\frac{1}{m}\\sum\\limits_{i=1}^{m}{({{h}}\\left( {{x}^{(i)}} \\right)-{{y}^{(i)}})x_{_{j}}^{(i)}}$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "def gradient(w, X, y):\n",
    "    w = np.matrix(w)\n",
    "    X = np.matrix(X)\n",
    "    y = np.matrix(y)\n",
    "\n",
    "    parameters = int(w.ravel().shape[1])\n",
    "    grad = np.zeros(parameters)\n",
    "\n",
    "    error = sigmoid(X * w.T) - y\n",
    "\n",
    "    for i in range(parameters):\n",
    "        term = np.multiply(error, X[:, i])\n",
    "        grad[i] = np.sum(term) / len(X)\n",
    "\n",
    "    return grad"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "注意，我们实际上没有在这个函数中执行梯度下降，我们仅仅在计算一个梯度步长。在练习中，一个称为“fminunc”的Octave函数是用来优化函数来计算成本和梯度参数。由于我们使用Python，我们可以用SciPy的“optimize”命名空间来做同样的事情。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "我们看看用我们的数据和初始参数为0的梯度下降法的结果。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([ -0.1       , -12.00921659, -11.26284221])"
      ]
     },
     "execution_count": 11,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "gradient(w, X, y)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "现在可以用SciPy's truncated newton（TNC）实现寻找最优参数。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(array([-25.16131872,   0.20623159,   0.20147149]), 36, 0)"
      ]
     },
     "execution_count": 12,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "import scipy.optimize as opt\n",
    "result = opt.fmin_tnc(func=cost, x0=w, fprime=gradient, args=(X, y))\n",
    "result"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "让我们看看在这个结论下代价函数计算结果是什么个样子~"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.20349770158947425"
      ]
     },
     "execution_count": 13,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "cost(result[0], X, y)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "接下来，我们需要编写一个函数，用我们所学的参数w来为数据集X输出预测。然后，我们可以使用这个函数来给我们的分类器的训练精度打分。\n",
    "逻辑回归模型的假设函数： \n",
    "\t\\\\[{{h}}\\left( x \\right)=\\frac{1}{1+{{e}^{-{{w }^{T}}X}}}\\\\] \n",
    "当${{h}}$大于等于0.5时，预测 y=1\n",
    "\n",
    "当${{h}}$小于0.5时，预测 y=0 。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [],
   "source": [
    "def predict(w, X):\n",
    "    probability = sigmoid(X * w.T)\n",
    "    return [1 if x >= 0.5 else 0 for x in probability]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "accuracy = 89%\n"
     ]
    }
   ],
   "source": [
    "w_min = np.matrix(result[0])\n",
    "predictions = predict(w_min, X)\n",
    "correct = [\n",
    "    1 if ((a == 1 and b == 1) or (a == 0 and b == 0)) else 0\n",
    "    for (a, b) in zip(predictions, y)\n",
    "]\n",
    "accuracy = (sum(map(int, correct)) % len(correct))\n",
    "print('accuracy = {0}%'.format(accuracy))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "我们的逻辑回归分类器预测正确，如果一个学生被录取或没有录取，达到89%的精确度。不坏！记住，这是训练集的准确性。我们没有保持住了设置或使用交叉验证得到的真实逼近，所以这个数字有可能高于其真实值（这个话题将在以后说明）。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 正则化逻辑回归"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "在训练的第二部分，我们将要通过加入正则项提升逻辑回归算法。如果你对正则化有点眼生，或者喜欢这一节的方程的背景，请参考在\"exercises\"文件夹中的\"ex2.pdf\"。简而言之，正则化是成本函数中的一个术语，它使算法更倾向于“更简单”的模型（在这种情况下，模型将更小的系数）。这个理论助于减少过拟合，提高模型的泛化能力。这样，我们开始吧。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "设想你是工厂的生产主管，你有一些芯片在两次测试中的测试结果。对于这两次测试，你想决定是否芯片要被接受或抛弃。为了帮助你做出艰难的决定，你拥有过去芯片的测试数据集，从其中你可以构建一个逻辑回归模型。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "和第一部分很像，从数据可视化开始吧！"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>Test 1</th>\n",
       "      <th>Test 2</th>\n",
       "      <th>Accepted</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>0.051267</td>\n",
       "      <td>0.69956</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>-0.092742</td>\n",
       "      <td>0.68494</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>-0.213710</td>\n",
       "      <td>0.69225</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>-0.375000</td>\n",
       "      <td>0.50219</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>-0.513250</td>\n",
       "      <td>0.46564</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "     Test 1   Test 2  Accepted\n",
       "0  0.051267  0.69956         1\n",
       "1 -0.092742  0.68494         1\n",
       "2 -0.213710  0.69225         1\n",
       "3 -0.375000  0.50219         1\n",
       "4 -0.513250  0.46564         1"
      ]
     },
     "execution_count": 16,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "path = 'ex2data2.txt'\n",
    "data2 = pd.read_csv(path, header=None, names=['Test 1', 'Test 2', 'Accepted'])\n",
    "data2.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAt8AAAHgCAYAAAB9zgEhAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAgAElEQVR4nOzdfXikZX33/883+IBmQ2UBW8yKrL9dK1pggWBFtkZUENcHwopmtSrarbTW3u6yh+0ut7/Wu7Ye2NW7afxpH+hKResNqzaJ9OhaHkSrK6IG70URkQSFugQLLj5MgqCQ7++Pcy5yZTIzmSQz1+P7dRw5JnNe12TOuXIl853z+p7n19xdAAAAADqvK+0OAAAAAGVB8A0AAAAkhOAbAAAASAjBNwAAAJAQgm8AAAAgIQTfAAAAQEIel3YHknT00Uf78ccfn3Y3AAAAUGA333zzj939mHrbShV8H3/88RofH0+7GwAAACgwM7u70TbSTgAAAICEEHwDAAAACSH4BgAAABJSqpxvAAAABL/61a908OBBPfTQQ2l3JbcOP/xwrVmzRo9//ONbfgzBNwAAQAkdPHhQPT09Ov7442VmaXcnd9xdhw4d0sGDB7V27dqWH0faCQAAQAk99NBDOuqoowi8l8nMdNRRRy35ygHBNwAAQEkReK/Mco4fwTcAAABSMzo6KjPT7bffnsjz/e3f/q0efPDBJT3mi1/8ol75yle25fkJvgEAALCoSkXas0fauTPcVirt+blXXnmlNm7cqKuuuqo9P3ARywm+24ngGwAAAE3t3y/19krbt0u7d4fb3t7QvhLT09P6yle+oo9+9KPzgu/du3frxBNP1Mknn6xdu3ZJkiYnJ/XSl75UJ598sk499VTdeeedkqQPfOADOv3003XSSSfpPe95jyTprrvu0rOf/WxdeOGFOumkk3TBBRfowQcf1Ic+9CFNTU3prLPO0llnnSVJuvbaa3XGGWfo1FNP1Wtf+1pNT09Lkv7jP/5Dz372s7Vx40aNjIys7IXGEHwDAACgoUpF2rQp3M7MhLaZmbn2aqy6LGNjYzr33HP1rGc9S6tXr9Y3v/lNfe5zn9PY2Ji+9rWv6ZZbbtGf/umfSpJ+93d/V+94xzt0yy236MYbb9Sxxx6ra6+9VhMTE/r617+uAwcO6Oabb9aXvvQlSdL3vvc9XXTRRfrWt76lI444Qn/3d3+nd77znXra056mL3zhC/rCF76gH//4x/qrv/orXX/99frmN7+pvr4+/c3f/I0eeughve1tb9O//du/6ctf/rJ+9KMfrfQwPobgGwAAAA3t3SvNztbfNjsbti/XlVdeqS1btkiStmzZoiuvvFLXX3+93vrWt+rJT36yJGn16tWqVCq65557dP7550sK62s/+clP1rXXXqtrr71Wp5xyik499VTdfvvtmpiYkCQ9/elP15lnnilJeuMb36j9dYbpb7rpJt12220688wztWHDBl1xxRW6++67dfvtt2vt2rVav369zExvfOMbl/8ia7DONwAAABqamJgb8a41MyNNTi7v5x46dEg33HCDbr31VpmZHn30UZmZXvOa1yxYRcTd6/4Md9cll1yiP/iDP5jXftdddy34GfVWJnF3nX322bryyivntR84cKBjK8Ew8g0AAICG1q+Xurvrb+vultatW97P/cxnPqM3v/nNuvvuu3XXXXfphz/8odauXavVq1fr8ssvf2xS5AMPPKAjjjhCa9as0djYmCTp4Ycf1oMPPqiXvexluvzyyx/L077nnnt03333SZL+67/+S1/96lclzU3qlKSenh5VqrNFn//85+srX/mKJqufIB588EHdcccdevazn60f/OAHj+WV1wbnK0HwDQAAgIYGB6WuBhFjV1fYvhxXXnnlY2kkkde85jWamprSq1/9avX19WnDhg364Ac/KEn6xCc+oQ996EM66aST9IIXvEA/+tGPdM455+gNb3iDzjjjDJ144om64IILHgusTzjhBF1xxRU66aST9MADD+jtb3+7JOmiiy7Sy1/+cp111lk65phj9LGPfUyvf/3rddJJJ+n5z3++br/9dh1++OG67LLL9IpXvEIbN27UM57xjOW9yDqs0TB+EfX19fn4+Hja3QAAAEjdd7/7XZ1wwgkt7bt/f5hcOTsbUk26u0PgvW+fVB1QzpS77rpLr3zlK3Xrrbd2/LnqHUczu9nd++rtT843AAAAmtq4UZqaCpMrJydDqsngoLRqVdo9yx+CbwCo5S6NjUkDA1J8wk2jdgAogVWrpK1b0+5Fa44//vhERr2Xg5xvAKg1NiZt3ixdfHEIuKVwe/HFob064QcAgKVi5BsAag0MSNu2ScPD4f7QUAi8h4dD+8BAuv0DAOQWwTcA1DILAbcUAu4oCN+2LbSTcgIAWCbSTgCgnngAHiHwzi93aXR0Lo1osXYA6BCCbwCoJ8rxjovngCNfyOMHMumwww7Thg0b9Fu/9Vt61atepZ/+9KeLPuYFL3jBsp5rbGxMt91225Ift6rNS7oQfANArSgoi3K8Z2fncsAJwPMpnscf/Q7J4wda08ErR0960pN04MAB3XrrrVq9erU+8pGPLPqYG2+8cVnPtdzgu90IvgGg1tjYXFAWpZoMDc0Fb4yS5k/t77Cra+HvGEB9CV05OuOMM3TPPfc8dv8DH/iATj/9dJ100kl6z3ve81h7fCS60T4f//jHddJJJ+nkk0/Wm970Jt144426+uqr9Sd/8ifasGGD7rzzTt15550699xzddppp+l3fud3dPvtt0uSfvCDH+iMM87Q6aefrj/7sz9ry2ubx91L83Xaaac5ACxqdtZ9ZCTcttKO/JiddQ9hQ/jid4kSu+2221rbcXbWfdu28DezbVv9+8vU3d3t7u6PPPKIX3DBBf65z33O3d2vueYaf9vb3uazs7P+6KOP+ite8Qr/z//8z3mPabTPrbfe6s961rP8/vvvd3f3Q4cOubv7hRde6J/+9Kcfe+4Xv/jFfscdd7i7+0033eRnnXWWu7u/6lWv8iuuuMLd3T/84Q8/9nyN1DuOksa9QTzKaicAUMtMOv/81tuRD43y+Bn5Bprr4ApQv/jFL7RhwwbdddddOu2003T22WdLkq699lpde+21OuWUUyRJ09PTmpiY0Atf+MLHHtton1tuuUUXXHCBjj76aEnS6tWrFzzv9PS0brzxRr32ta99rO3hhx+WJH3lK1/Rv/7rv0qS3vSmN2nnzp3Lfn31kHYCACg+8viBlenQClBRzvfdd9+tX/7yl4/lfLu7LrnkEh04cEAHDhzQ5OSkttaU12y0j7vLFunX7OysnvKUpzz22AMHDui73/1u7OV27gM5wTcAoPjI4wdWpsMrQP3ar/2aPvShD+mDH/ygfvWrX+llL3uZLr/8ck1PT0uS7rnnHt13333zHtNon5e85CX61Kc+pUOHDkmSHnjgAUlST0+PKpWKJOmII47Q2rVr9elPf7r68ly33HKLJOnMM8/UVVddJUn65Cc/2ZbXF0fwDQAovoEBaWRk/khdFICPjLDaCdBMQleOTjnlFJ188sm66qqrdM455+gNb3iDzjjjDJ144om64IILHguco1HpRvs897nP1bvf/W719/fr5JNP1o4dOyRJW7Zs0Qc+8AGdcsopuvPOO/XJT35SH/3oR3XyySfruc99rj772c9KkoaHh/WRj3xEp59+un72s5+15bXFmZfoUltfX5+Pj4+n3Q0AAIDUffe739UJJ5yw+I6jo2FVk/iVo3hAPjKS2HyYQ4cO6dRTT9Xdd9+dyPO1ot5xNLOb3b2v3v5MuAQAAEBj0ZWjgYGFV476+xO7cjQ1NaUXvehFete73pXI83UKwTcAAAAay8gKUE972tN0xx13JPZ8nULONwAAAJAQgm8AAICSKtPcv05YzvEj+AYAACihww8/XIcOHSIAXyZ316FDh3T44Ycv6XHkfAMAAJTQmjVrdPDgQd1///1pdyW3Dj/8cK1Zs2ZJjyH4BgAAKKHHP/7xWrt2bdrdKB3STgAAAICEpBp8m9nlZnafmd3aYLuZ2YfMbNLMvmVmp8a2XWhmE9WvC5PrNUrBPRQVqM2Da9QOtIpzCwBKLe2R749JOrfJ9pdLWl/9ukjS30uSma2W9B5Jvy3peZLeY2ZHdrSnKJexsVDNK142N6rmtXlz2A4sB+cWAJRaqsG3u39J0gNNdjlP0sc9uEnSU8zsWEkvk3Sduz/g7j+RdJ2aB/HA0gwMhDK6w8NzQVJURnfbtsSqeaFD0hx95twCgFLL+oTLXkk/jN0/WG1r1A60R1Q2VwpB0fBw+H7bttAelddFPkWjz/HfZzwIHhnpXNU2zi0AKLW0004WU+9dyJu0L/wBZheZ2biZjbOUDpYkHiRFCI6KIe3RZ84tACitrAffByU9PXZ/jaSpJu0LuPtl7t7n7n3HHHNMxzqKAooCsrh4ni7yKwp+owC8q2su8E4iCObcAoDSynrwfbWkN1dXPXm+pJ+5+72SrpF0jpkdWZ1oeU61DWiP2pHQ2dmFI6XIt7RGnzm3AKDUUs35NrMrJb1I0tFmdlBhBZPHS5K7/4OkfZI2SZqU9KCkt1a3PWBmfynpG9Uf9V53bzZxE1iasbGFI6HxPN3+/s7lBCMZjUafOx2Ac24BQKmZl2iUpa+vz8fHx9PuBvLAPQRJAwPzA7FG7SVSqUh790oTE9L69dLgoNTTk3avlqh29HloaOH9Tv1+ObcAoPDM7GZ376u7jeAbQKv275c2bQqZEjMzUnd3SJfet0/auDHt3i3B6Gh6q50AAAqP4LuK4BtYvkpF6u0Nt7V6eqSpKWnVquT7tSyMPgMAOqhZ8J31CZcAMmLv3jDiXc/sbNieG2ZhZLs2wG7UDgBAmxB8A2jJxERINalnZkaanEy2PwAA5BHBN4CWrF8fcrzr6e6W1q1Ltj8AAOQRwTeAlgwOhsmV9XR1he1AXe5hkmvtHKNG7QBQYATfAFrS0xNWNenpmRsB7+6ea8/NZEskb2wsrC4TLyIUrS6zeXPYDgAlkWqRHQD5snFjWNVk796Q471uXRjxJvBGUwMDc1U8pYXrqg8MpNs/AEgQSw0CWcRSeGhVXs6V+DrqkU4XNAKAlLDUIJA3XKZHq/JyrpiFQDuOwBtACRF8A1kUv0wfBVVcpkc9eTlXon7FxT8wAEBJEHyjnLK++kI0ShgFVV1dc8EUo4WIy8O5UvuBYHZ24QcGACgJcr5RTqOj4ZJ8PECJBwgjI6HSYdrc56/vNzubjWAK2ZPlcyUvf28A0CbkfAO18nCpnsv0aFXWz5WBgRBgx0fioxH7kZFs/L0BQEIIvlFOWb9Uz2V6tCoP54pZGNmu/btq1A7EZT1NEFgigm+UV5ZXXxgbW/hhIP5hISsrWCB9nCsourys6AO0iOAb5ZXlS/VcpkerOFdQdHlIEwSWgAmXKKfaf961FfeyMgIOAKBIE3Kn2YRLgm+UE6svAEC+ZHlFH6AGq50AtbhUDwD5keU0QWCJCL5RTqy+AAD5kIcVfYAleFzaHQAAAGio0Yo+Umjv7ydNELlC8A0AALIrShMcGFiYJtjfT5ogcofgGwAAZFeUDthqO5Bx5HwDAAAACSH4BgAAABJC8A0AAAAkhOAbAAAASAjBNwAASXEPFXZr16Zu1A6gcAi+AaANKhVpzx5p585wW6mk3SNk0tiYtHnz/OIwURGZzZvDdgCFxlKDALBC+/dLmzaFwnszM1J3t7Rjh7Rvn7RxY9q9Q6YMDMxVZ5TCWtXx6o2sWQ0UnnmJLnH19fX5+Ph42t0AUCCVitTbW3+ku6dHmpqSVq1Kvl/IsHi59Ei8eiOA3DOzm929r9420k4A5E6WUjz27g0j3vXMzobtwDzx8ugRAm+gNAi+AeTK/v1hpHn7dmn37nDb2xva0zAxEVJN6pmZkSYnk+0PciAa+Y6L54ADKDSCbwC5UamE3OpKZS7gnZmZa5+eTr5P69eHHO96uruldeuS7Q8yLp5ysm1buDwS5YATgAOlQPANIDeymOIxOCh1NfhP2tUVtgOPGRubC7yjVJOhobkAnNVOgMIj+AaQG1lM8ejpCaua9PTMjYB3d8+1M9kS8wwMSCMj83O8owB8ZITVToASIPgG8JgsTWSsJ6spHhs3hlVNhoelXbvC7dRUDpYZzHLBlyz3bSXMpPPPXzi5slE7gMIh+AYgKXsTGevJcorHqlXS1q3SpZeG21yMeGe54EuW+wYAK0DwDSCTExnrIcWjzeIFX6IgNysFX7LcNwBYASpcAmhpIuPWrcn2qZEoxWPv3pDjvW5dGPEm8F6G+HrTw8NzRV+yUPAly30DgBWgwiUA7dwZUk0a2bUrpFOgoNzn5/PMzmYnuM1y3/LOPaTvDAzMP6aN2gG0jAqXAJrK6kRGJCDLBV+y3LciIK8eSEWqwbeZnWtm3zOzSTPbVWf7kJkdqH7dYWY/jW17NLbt6mR7DhRLlicyooOyXPAly30rCvLqgVSklvNtZodJ+oiksyUdlPQNM7va3W+L9nH3i2P7/w9Jp8R+xC/cfUNS/QWKLJqwuGlTiHFmZsKId1cXExkLrVHBFym09/eH5e/oWzGRVw+kIrWcbzM7Q9L/cveXVe9fIknuXjez1MxulPQed7+uen/a3ZcUEpDzDTQ3Pc1ExlLJcs5vlvtWNOTVA23XLOc7zdVOeiX9MHb/oKTfrrejmT1D0lpJN8SaDzezcUmPSHq/u5OcBqxQtFY1SiIq7NJqe5Ky3LciaZRXz8g30DFp5nzX+6tuNAy/RdJn3P3RWNtx1U8Ub5D0t2b2/9R9ErOLzGzczMbvv//+lfUYAICiyFJefVErmgJ1pBl8H5T09Nj9NZKmGuy7RdKV8QZ3n6refl/SFzU/Hzy+32Xu3ufufcccc8xK+wwAQDE0yquPAvAkVzth5RWUSJrB9zckrTeztWb2BIUAe8GqJWb2m5KOlPTVWNuRZvbE6vdHSzpT0m21jwUAAA0MDEgjI/NTTKIAfGQk2dVOWHkFJZJazre7P2JmfyzpGkmHSbrc3b9jZu+VNO7uUSD+eklX+fyZoSdI+kczm1X4APH++CopAABgEVnKq2flFZQIFS4BAEA2sPIKCoIKlyi2vE7UyWu/AaATqGiKkiD4Rv7ldaJOXvudgkpF2rNH2rkz3FYqafcoWzg+yL0srbwCdJq7l+brtNNOcxTQ7Kz7tm3uUritdz+L8trvhH35y+49Pe7d3eHQdHeH+1/+cto9ywaODwphZGTh/774/8SRkXT7ByyRwvzFuvEoOd8ohvioSSQPE3Xy2u8VqlRCJc2JCWn9+lBJs6en/n69vfVHcnt6pKmpclfg5Ph0VqvnKdqAiqYomGY53wTfKI68TtTJa7+Xaf9+adOm8DJnZqTu7vDy9+2TNm6cv++ePdL27WG/Wt3d4TNLmStycnw6ZynnKQDUYsIlFpf3yX95naiT134vU6USAppKZS5gnJmZa5+enr//xET9wDJ63ORk+/uXp9zppI9PWSz1PAWApSD4RpDnyX95naiT136vwN694WXWMzsbtsetXx9GHOvp7pbWrWtf3/bvDykc27dLu3eH297e0J5VSR6fMlnqeQoAS0HwjSDP1cWyVCJ5KfLa7xVY6kjt4OD8jJy4rq6wvR3yOtKZ1PEpG64oAOik1CpcImPyXF0sKpEcn5ATvZ7+/ux+cMhrv1cgGqltlKNcO1Lb0xNybBvl3rZrMmErI51ZzJ1O6viUzVLP06QxERTINyZcYr6STf5Dspa7Osf0dAg2JidD4DM42N7AcufOkGrSyK5d0qWXtu/52q3Tx6dssryKDBNBC45VXwqj2YRLRr4xp9Hkv6yPfCM3ljtSu2pVZ0eesz7SuZhOH5+yyeoVhXh6VCQ6ZzdtYmnJQojmX8WvOsfTQEdGpPPPT7uXWCGCbwS1Od5DQ/PXnyYAR5ts3BiChCyN1A4OSjt21N9G7nQ5ZfE8zWt6FJYgPv9Kmv9enPX5V2gZwTeCRpP/pNDe38+nbbTNUkdqO53jmtWRTqQra1cUmAhaAnmef4WWkfONgDwzZFSSOa7kTiPLKKpUIsy/yj0qXFYRfAP5kuWJb0DS+HtoQREGkuJpoBFGvnOHCpfACuSt6mGRUOwEmBOlR/X0zBVX6u6eay994C3lu2CcVMria2VEzjfQRL2Uhx07WNYrKeS4AvNlcSJopuR9wiLzr0qB4BtogGW90pf3JQCBTsjaRNBMyfuExRIWXysj0k6ABkh5SB/l09ER7tLo6MJL+I3akS/xADySh8BbCn08//yFfW3Ujlwi+AYaIOUhfeS4oiPynheM5hoVjONDFTKCtBOgAVIesoEcV7Rd3vOC0RgF45ADLDUINMCyXkCBsZxbMY2OUp4dmcA631UE31iqJAu8AEgYhUyKpwjrfKMQWOcbmZbldbSjlIfhYWnXrnA7NUXgDeQeecHp6tSkVyYsIgfI+Uaq8rCONst6AQVDXnD6okmvpIeghAi+kRrW0QaQCgqZpI9Jrygxgm+kppV1tBlxBtB2LRYyqVTC/6GJibD60eBgmGyNNsh7MRxgBZhwidTs3Cnt3t14+65d0qWXJtcfAIsrS0DKZOuEMOkVBcWES2RStI52PayjDWTP/v1h+c3t28MH5+3bw/39+9PuWXvFU+KiVLiZmbn26el0+1cYTHpFSRF8IzWUDgfyo0wBaSspcVih2kmvs7NzOeAE4Cg4gm+khtLhQH6UKSCdmKhf2VYK7ZOTyfankBpNeo0C8LGxtHsIdAwTLpEqSocD+VCmgDRKiav3ekmJa5MWJ70CRUTwjdSxjjayoiyTCZejTAHp4GCoN1APKXFtEhW9abUdKBBWOwFQKMsNoFndorlKJUyurFeBtqeneOvycz4AWIlmq50QfAMojOUGTGULLJerbAHp9DQpcQCWp1nwTdpJmbiHSSzxHLtm7UCOrKRiKgWfWtPqHI2ipO+QEgegEwi+y2RsTNq8ef7s8vhyTyMj5Noht1YSQJdpMuFKLRaQ1hsd37GjuKPjALBULDVYJgMDC9dRja+zyuxy5NhKAmgKPrVHmdYCB4DlIvguk9p1VLu6Fq6zCuTUSgJoCj61RyfWAq9UpD17pJ07w229vHwAyBOC77KJAvA4Am8UwEoCaAo+tUe703fKUs4eQLkQfJdNlGoSRylfFMBKA+hoMuHwsLRrV7idmiJPeSnamb5DCguAomLCZZnU5ngPDc3dlxgBR+6ttGIqq1usTDuL07ACDYCiIvguk7GxhTneUQrK8HAo6ctqJ5lUlKXbkkAAnZ7oKkOjtcCXkr7DCjTIDZbxxRKlmnZiZuea2ffMbNLMdtXZ/hYzu9/MDlS/fj+27UIzm6h+XZhsz3NqYCAsJxgf4Y4C8JERVjvJKPJekSftSt9hBRrkRrSMbzyFM7rSvHlz2A7EpFbh0swOk3SHpLMlHZT0DUmvd/fbYvu8RVKfu/9xzWNXSxqX1CfJJd0s6TR3/0mz56TCZcYxerAAlRdRVpz7yI1mKZ2sJlZazSpcpjny/TxJk+7+fXf/paSrJJ3X4mNfJuk6d3+gGnBfJ+ncDvUTSWH0YIFOLN0G5AEr0GQbS0DGsIwvlijNnO9eST+M3T8o6bfr7PcaM3uhwij5xe7+wwaP7e1UR5GQeBEgaeHoQQnTYsh7RZmtdAItOoMqpnVEAXj0/iUReKOhNIPvemdkbQ7Mv0m60t0fNrM/lHSFpBe3+NjwJGYXSbpIko477rjl9xadVzsBNPonVuLRgyjvtV4ATt4ryoAJtNkSXwIyEv1/2rSpxOlAjZbxLel7F5pLM+3koKSnx+6vkTQV38HdD7n7w9W7/yTptFYfG/sZl7l7n7v3HXPMMW3pODqIIkDzpF15kUvLAOJIhaujNud7dnbuKi51NFBHmsH3NyStN7O1ZvYESVskXR3fwcyOjd19taTvVr+/RtI5ZnakmR0p6ZxqG/KOIkDzpJn3yiorAGqRCldHo2V8owC8hPOV0FxqaSfu/oiZ/bFC0HyYpMvd/Ttm9l5J4+5+taR3mtmrJT0i6QFJb6k+9gEz+0uFAF6S3uvuDyT+ItBeCRYBytO62WnkvXJpGUA9pMLVES3jG1+RKwrA+/tLOV8JzaW21GAaWGow40ZHw6om8dGDeEA+MtKWIkD1JgtFRUBKO1moxp49YaS70Rvs8DB5uEAZsQQk0JqsLjUIzJdAEaD4iG4UWM7MzLVPT6/4KQqBS8sA6mEJSGDlKC+P7DCrP7LdqH0ZWpksxIgul5YBNMYSkMDKEHyjVBjRbc3gYFi3t54kVlkBkG0sAQksH8E3SoUR3dZEl5Ab5cYzwgWg6PI0MR/5woRLlAqThZZmeppLywDKh4n5WKlmEy4Z+UapMKK7NFxaBlA2LLWKTiP4RukwWQgA0AgT89FpBN8oJUZ0AQD1MDEfncY63wAAAFXRxPx6mJiPdiD4BgAAqBocDPOA6mGpVbQDwTcAAEAVVTzRaeR8AwAAxDAxH51E8A0AKaCAB5BtTMxHpxB8A0DC6hXw2LGDAh4AUAbkfANAguIFPKLlzGZm5tqnp9PtHwCgswi+ASBBrRTwAAAUF8E3ACSIAh4AUG4E3wCQIAp4AEC5EXwDQIIo4AEA5UbwDQAJooAHgExzl0ZHw20r7Vgygm8ASFhUwGN4WNq1K9xOTbHMIIAMGBuTNm+WLr54LtB2D/c3bw7bsSKs8w0gFWUvMkMBDwCZNDAgbdsWRgUkaWgoBN7Dw6F9YCDd/hWAeYkuH/T19fn4+Hja3QBKr16Rma4uiswAQCZEI91RAC6FwHtoSDJLr185YmY3u3tf3W0E3wCSVKlIvb3htlZPT0i/IO8ZAFLmPn92+OwsgfcSNAu+yfkGkCiKzABAxkUj33HxHHCsCME3ElGpSHv2SDt3htt6o54oB4rMAECGxVNOtm0LoyJRDjgBeFsw4RIdVy+/d8cO8nvLKioyUy8Ap8gMAKRsbGwu8I5yvIeGwrbhYam/Xzr//HT7mHPkfKOjyO9FLc4JAMgw9xCADwzMz/Fu1I66yPlGasjvRS2KzABAhpmFke3aALtRO5aMtOxQ8W8AACAASURBVBN0FPm9qCcqMrN3bzgH1q0L63wTeAMAio7gu9NKfvmG/N7OKEKBGorMAADKiLSTTit5mdbBwfnLhMZ1dYXtWJr9+0PO9Pbt0u7d4ba3N7QDAIBsI/jutHiZ1igAL1GZVvJ726tSCSvHVCpzVxNmZubap6fT7R8AAGiOtJNOq12iJyrVWqIyreT3tk8rE1hJ5QAAILtaGvk2szVmdlb1+yeaWXdnu1Uw8QA80q7A210aHV246H2j9pRE+b2XXhpuCbyXhwmsAADk26LBt5n9nqSrJe2pNj1D0mc72anC6WSZ1pLnlJdNNIG1HiawAgCQfa2MfL9T0vMl/VyS3P0OSU/tZKcKpdNlWkueU142TGAFACDfWgm+H3L3X0Z3zOwwScVPVG6XRmVao4B5pSPTtT+vq2vh86EwmMAKAAWWk1RSrMyi5eXN7H9L+m9Jb5X0R5LeIWnC3S/pfPfaK5Xy8kmt8+0+f0h0dpbAu8Cmp5nACgCFMzoaUkbjA2jxK9ojI6HKJDKvWXn5VlY7+VNJF0m6XdI2SddI+sf2da/gonKsrbYvR6Occka+C4sCNQCyrgjFwBIXTyWVwvs4qaSF0zT4rqaYXO7uF0r6+2S6hCWpzfGO/6FKBOAAgMTt3x9qD8zOhpWYurulHTtCetzGjWn3LsNYnrgUWkk7uVbSK9z9V8l0qXNSSTvpNC5RAQAypFIJVXcrlYXbenpC3QfS5BZBKmnuNUs7aWXC5fclfdnMLjGzd0Zf7e0ilm1gIATY8U/E0SfnkREuUQEAEtVKMTA00cnliZEJrQTf90u6TtKTJR0T+1oxMzvXzL5nZpNmtqvO9h1mdpuZfcvMPm9mz4hte9TMDlS/rm5Hf3Ipyh2v/UTcqB0AgA6iGNgKdHp5YmTCohMu3f3PJMnMnlS9/4t2PHE1n/wjks6WdFDSN8zsane/Lbbb/5XU5+4PmtnbJe2WFK1k/At339COvgAAgPaIioHVC8ApBraIRssTS6G9v59U0gJopcLlc8zsG5ImJE2a2dfM7IQ2PPfzJE26+/er64hfJem8+A7u/gV3f7B69yZJa9rwvAAAoEMoBrYCpJKWQitpJ5dJ+p/uvsbdeyW9W9I/teG5eyX9MHb/YLWtka2SPhe7f7iZjZvZTWbG2QgAQAZQDGwFSCUthVbW+e5x9+uiO+5+fbXwzkrVO4PqJjOZ2Rsl9UnqjzUf5+5TZvZMSTeY2bfd/c46j71IYZ1yHXfccSvvNQAAaGrjxrCqCcXAgIVaCb7vMrNLJH2iev+Nku5uw3MflPT02P01kqZqdzKzlyqMtve7+8NRu7tPVW+/b2ZflHSKpAXBt7tfpjB6r76+PmYqAACQAIqBAfW1knbyewpB8r7q1xqFUvMr9Q1J681srZk9QdIWSfNWLTGzUxSqab7a3e+LtR9pZk+sfn+0pDMlxSdqAgAAAJnTymonhyT9Ubuf2N0fMbM/VihXH1XS/I6ZvVfSuLtfLekDklZJ+rSFPKf/cvdXSzpB0j+a2azCB4j316ySAqSKssoAAKCeVipc/oekLe7+0+r9IyX9i7u/IoH+tVUhK1yiZUkFxPXKKnd1UVYZAICyaFbhspWc71+PAm9JcvefmNnT2tY7IAH1AuIdO9ofEFcq4XniZZWjtW43baKsMgAAZddKzvesmT22vraZsWQIciUeEEeB8MzMXPv0dPuei7LKAACgmVaC7z+X9BUz+2cz+2dJX5L0PzvbLaB9kgyIKasMAACaaWXC5b+b2fMknVFt2hlfeQTIuiQDYsoqAwCAZhqOfJvZ083sCEly9/+W9ICk35G0xcwen1D/gBWLAuJ62h0QU1YZAAA00yzt5NOSjpAkMztZ0qik+yQ9T9JHOt81oD2SDIgpqwwAAJpplnbyZHc/WP3+jQrrcP+1mXVJuqXzXQPaIwp8Gy3/1+6AmLLKAACgkWbBt8W+f7FCiXe5+6yZUaYduZJ0QExZ5cVRiAgAUEbNgu//NLP/I+leSUdJukGSzOw3JP0qgb4BbUVAnB1JrbsOAEDWNMv5fqekfZJ+JOl33P2X1fanSfqzTncMQDElue46AABZ03Dk291nJf1LnfZvdrRHAAqtlXXXuUIBQCI9DcXUSnl5AGgbChEBaAXpaSiqVipcAkDbJLnuOoB8Ij0NRUbwDSBRFCICsJhW0tOAvGpW4bLHzP7SzP7ZzF5Xs+3/63zXABQRhYgALIb0NBRZs5zvyyXdLenfJf2emb1G0hvd/VeSzkyicwCKiUJEAJqJ0tPqBeCkpyHvzL1+vRwzO+DuG2L33yPppZJeLenz7n5qMl1sn76+Ph8fH0+7GwAAoIlKRertDbe1enrCh3c+rCPLzOxmd++rt61Zzvfh1VLykiR3/wtJH5P0ZUmr29pDAACAKtLTUGTN0k7+XdJLJF0XNbj7R83sR5I+3OmOAQCA8iI9DUXVMO2kiEg7AQAAQKctN+0EAAAAQBsRfAMAAAAJWTT4NrMFeeH12gAAAIDEuUujo+G2lfaUtTLy/fUW2wAAAIBkjY1JmzdLF188F2i7h/ubN4ftGdJwBNvMnirpWElPMrMTJVl10xGSnpxA3wAAAIDmBgakbduk4eFwf2goBN7Dw6F9YCDd/tVolj7yCkm/J2mNpI9oLviuSPqzDvcLAAAAWeYeRpUHBiSzxds7xSwE3FIIuKMgfNu20J5EH5Zg0aUGzex17v6phPrTUSw1CAAA0CajoyGtIx7kRukew8PSyIh0/vnJ9cdd6oplVM/OphZ4r3Spwaea2RHVH/QPZvZ1M3tJW3sIAACAfImne0T51mmle0TPHRfPAc+QVoLvi9z952Z2jkIKytsl7e5stwAAAJBpUbpHFIB3dc0F3kmme9QG/bOzCz8UZEgrwXfU45dL+md3v7nFxwEAAKDI4vnWkaTzrMfGFgb98Q8FGVvtpJUg+hYz2yfpVZI+Z2arNBeQAwAAoKyykO4xMBDyy+NBfxSAj4xkbrWTVoLvt0r6X5Ke5+4PSjpc0tZOdgoAAAAZl5V0D7MwsbN2tL1Re8oWrVTp7o+a2TMlnS3pfZKeJNJOAAAAyq1RuocU2vv7k13tJCdaWWrww5IeL+mF7n6Cma2WdI27n55EB9uJpQYBAADaJCvrfGfQSpcafIG7/4GkhyTJ3R+Q9IQ29g9Z4B7W66z9MNaoHQAAlFvO0j2yopXg+1dm1qXqJEszO0rSbEd7heSNjYWF8uM5WlEu1+bNmZspXESVirRnj7RzZ7itVNLuEQAAaLeGOd9m9jh3f0ShtPy/SjrGzP5C0usk/UVC/UNS4gvlSyFnK62F8kto/35p06YwV2VmRurulnbskPbtkzZuTLt3AACgXRrmfJvZN9391Or3z5X0Ukkm6Xp3vzW5LrYPOd+LiM9ajiS9UH4JVSpSb2/9ke6eHmlqSlq1qvN92LtXmpiQ1q+XBgfDcwMAgKVrlvPdLPj+v+5+Skd7ljCC7xa4hwpVkdlZAu8O27NH2r49jHjX6u4On4W2dnBxz3qj7l1djLoDALBczYLvZksNHmNmOxptdPe/WXHPkC2NFspn5LujJibqB95SaJ+c7NxzVyoh8I6Pukd92bQpmVF3AADKpNmEy8MkrZLU0+ALRZKVhfJLaP36MNpcT3e3tG5d5557797wq65ndjZsB5APTNpG27ESWkc0G/m+193fm1hPkK4MLZRftvzjwcEwubKerq6wvVPSHHUH0D5M2kZHRCuhxWOD+GDdyAhFdJahWfBNnkGZDAyEP6L4gvhRAN7fn9hqJ2V8A+npCa+vUd51J9M+olH3RvnmnRx1B9AepI+hY1gJrSOaTbhcXS2o07knNztX0rBCissed39/zfYnSvq4pNMkHZI06O53VbddImmrpEclvdPdr1ns+ZhwmW1ZWPUjTdPTYcR/cjIEvYODyaxyUuZjDhRB2pO2UXCshLYsy6pwmUDgfZjCGuIvl/QcSa83s+fU7LZV0k/cfZ2kIUl/XX3scyRtkfRcSedK+rvqz0OOlT3/eNWq8AZ56aXhNomgNxp17+mZyzvv7p5rJ/AGso/0MXRUPA01QuC9Iq1UuOyU50madPfvu/svJV0l6byafc6TdEX1+89IeomZWbX9Knd/2N1/IGmy+vOQY7yBpGPjxjDCPTws7doVbqemipvmAxRNmpO2UQKNVkJjsuWyNcv57rReST+M3T8o6bcb7ePuj5jZzyQdVW2/qeaxvZ3rKpJA/nF6olH3TivbZFogCWlO2kbB1a6EFs/5lhgBX6Y0R77r/bZqP0Y12qeVx4YfYHaRmY2b2fj999+/xC4iSYOD8+v7xPEGkn/794f88u3bpd27w21vb2gHsHykj6FjGq2EFk3CHBtLu4e5lObI90FJT4/dXyNpqsE+B83scZJ+TdIDLT5WkuTul0m6TAoTLtvSc3REmqt+oLNYjQHorCh9LOlJ2yi4jKyEVjRpBt/fkLTezNZKukdhAuUbava5WtKFkr4q6QJJN7i7m9nVkv6Pmf2NpKdJWi/p64n1HB3DG0gxtTKZltUYgJVJKn0sMe5hZDUe+DVrR/uZ1V/Hu1E7WpJa8F3N4f5jSdcoLDV4ubt/x8zeK2nc3a+W9FFJnzCzSYUR7y3Vx37HzD4l6TZJj0h6h7s/msoLQdsV7g0ETKYFsHQUeEFBpTnyLXffJ2lfTdufx75/SNJrGzz2fZLe19EOAmgLJtMCWDIKvKCgGhbZKSKK7ADpoJgPgGWhwAtyallFdgCgXViNAcCyUOAFBZRq2gmA8mAyLYAla1TghQAcOcbIN4A57tLo6MLKZY3alyiaTHvppeGWwBtAQ7UFXmZn53LAqbCIHCP4BjAnWl0g/sYWvQFu3kxBBQDJocALCorgG8nr8OgqViC+ukAUgLO6AIA0RAVe4ikmUQAeFX4BcojgG8ljdDW7akeWuroWjjwBQBKiQi61/3catQM5QfCN5DG6mm2sLgAAQMcQfJdB1tI8GF3NtkarC5AOBADZk7X3eCyK4LsMspjmwehqNrG6AADkSxbf49EUwXcZZDHNg9HVbGJ1AQDIlyy+x6MpysuXRZZK9Nb+YxgaWnifEfB0uIcAe2Bg/u+gUTsAIH1Zeo+HpObl5Qm+y8Q95FdHZmfT+aMcHQ2XwuL/GOL/OEZGwkx2AADQmqy8x0NS8+CbtJOyyFKaB2u3AgDQPll6j8eiCL7LIGuT6Fi7FQCA9sjaezwW9bi0O4AENJpEJ4X2/n7SPAAAyCPe43OHnO8yYBIdAADFxHt8JjHhsqq0wTcAAAASw4RLoB2oIgYAAFaI4BtoFVXEAADACjHhEmhVvIqYtLA4EEskAkDHVSrS3r3SxIS0fr00OCj19KTdK6B15HwDS0EVMQBIzf790qZNYTW9mRmpuzvUldm3T9q4Me3eAXOYcFlF8I22oIoYACSuUpF6e8NtrZ4eaWpKWrUq+X4B9TDhEmgXqogBQCr27g1jHfXMzobtQB4QfAOtoooYcqBSkfbskXbuDLf1RgmBPJqYCKkm9czMSJOTyfYHWC4mXAKtoooYMq5ePuyOHeTDohjWrw/ndL0AvLtbWrcu+T4By0HON9Aqqoghw8iHRdFxjiNPyPkG2sEsjGzXBtiN2oEEkQ+LouvpCVdxenrCSLcUbqN2Am/kBWknAFAA5MOiDDZuDCPce/eGc3rdurDON4E38oTgGwCWIKsFPsiHRTNZPW+XY9UqaevWZJ+zSMcP6SPnGygx3lCWJssFPsiHRSNZPm/zgOOH5aDIThXBNzCHN5SlyUNwy+8UtfJw3mYZxw/LxYRLAPNUKiFIq1Tm0hRmZubap6fT7V8W5WFCY5QPOzws7doVbqemCLzLLA/nbZZx/NAJBN9ACfGGsnR5mdAY5cNeemm4ZVSu3PJy3mZV4Y6fuzQ6urAoXKN2dATBN1BChXtDSUA0obEeJjQiqzhvV6Zwx29sTNq8eX5V5qh68+bNYTs6juAbS8cn59zr9BtKEUucDw6G/Ol6urrCdiBrOG9XpnDHb2AgVGkeHp4LwC++eK5688BA2j0sBYJvLB2fnHOvk28o+/eHCUrbt0u7d4fb3t7QnmcU+EAecd6uTOGOn5k0NDQXgHd1zQXeQ0MUi0sIq51g6Wo/KQ8NLbzPH3DmdWJljDKsDDA9TYEP5A/n7coU7vi5zx+BmZ3lfbvNWGqwiuC7jeIBeITAO3fa/YayZ08Y6W5U6GV4OPniGEAt1rdHqfH+nYhmwTcVLrE80aWr+B8vf7i50+5KcUzkRNbVu+KzYwdroaMkml25lngfTwg531ie6A84Lp4DjlIq3MoAKBTWt0fpjY0tTBGN54AzZysRBN9YutpPzrOzC2dPo5QKtzIACoX17fOtiKsoJW5gQBoZmT/CHQXgIyOsdpIQ0k6wdI0+OUuhvb9fOv/8dPuIVEQrADSayJnrCUrIPdKi8ot0oTYxq//+3KgdHZFK8G1mqyXtlXS8pLskvc7df1KzzwZJfy/pCEmPSnqfu++tbvuYpH5JP6vu/hZ3P5BE36G5T84DAws/Off388m55KIS54VaGQCFEKVFNZoQTFpUNsXThSLR73DTpmKsooRySSvtZJekz7v7ekmfr96v9aCkN7v7cyWdK+lvzewpse1/4u4bql8E3lJyxW+iT8i1kzIataM1BSpeRIlzZBFpUflEuhCKJq3g+zxJV1S/v0LSgqFSd7/D3Seq309Juk/SMYn1MI8ofpNv/P6AjipcwZSSIF0IRZNWzvevu/u9kuTu95rZU5vtbGbPk/QESXfGmt9nZn+u6si5uz/c4LEXSbpIko477rh29D274mVjpYXFb0gHyTZ+f0DHkRaVP6QLoWg6VmTHzK6X9Bt1Nr1b0hXu/pTYvj9x9yMb/JxjJX1R0oXuflOs7UcKAfllku509/cu1qdSFNlh8fxMW7S4B7+/FaOAClAsZaici+LJXIVLM/uepBdVR72PlfRFd//NOvsdoRB4X+run27ws14k6V3u/srFnrcUwbdE2diMarmcO7+/ZWv5GAPIFf62kTfNgu+0cr6vlnRh9fsLJX22dgcze4KkUUkfrw28qwG7zMwU8sVv7Whv84TiN5nUcnEPfn/LRgEVoLiidKHhYWnXrnA7NUXgjXxKK/h+v6SzzWxC0tnV+zKzPjPbU93ndZJeKOktZnag+rWhuu2TZvZtSd+WdLSkv0q2+xlF8ZvMamm2Pr+/FWFFBKDYWEUJRZHKhEt3PyTpJXXaxyX9fvX7f5H0Lw0e/+KOdjCvKH6TWS3N1uf3tyKsiAAAyAMqXBYJxW8yq6XZ+vz+VoQVEQAAeZDKhMu0lGbCJZbPPYxAxwPgZu0tYrZ+53GMAQBZkcUJl0A2dajQDcU9Oo9jDADIA9JOgLgOFrqhuEfncYwBAFlH2glQi0I3AABgBTJXZCctBN9oGYVuAADAMpHzDSwFhW4AAECHEHwDcRS6AQAAHcSESyCOQjcAAKCDCL6BOArdAACADiL4BuLM6o9sN2oHAABYAnK+AQAAgIQQfAMAAAAJIfgGAAAAEkLwDQAAACSE4BtAfrlLo6ML119v1A4Aecf/vdwj+AaQX2Nj0ubN8wsgRYWSNm8O2wGgSPi/l3ssNQggvwYG5iqQSmE99niFUtZlB1A0/N/LPfMSXZ7o6+vz8fHxtLsBoJ2iEZ/ojUiaX6EUAIqG/3uZZ2Y3u3tf3W0E3wByz13qimXRzc7yBgSg2Pi/l2nNgm9yvgHkWzQCFBfPhQSAouH/Xq4RfAPIr/il123bwshPlAvZ4I2oUpH27JF27gy3lUoK/QaA5VrG/z1kCxMuAeTX2NjcG1CU6zg0FLYND0v9/dL55z+2+/790qZN4b1qZkbq7pZ27JD27ZM2bkzpNTRQqUh790oTE9L69dLgoNTTk3avAKRuif/3kD3kfAPIL/fwRjQwMD/XsU57pSL19tYf6e7pkaampFWrEur3Iup9SOjqyuaHBAAJW8L/PaSHnG8gSyiQ0D5mYYSn9o2mTvvevSGYrWd2NmzPgkolBN6VSgi8pXAbtU9Pp9u/NJAqBMQs4f8esongG+krWzBKgYRUTEzMBbO1Zmakyclk+9NIXj4kJGX//nDFYvt2affucNvbG9oRU7b/o0COEXwjfWULRuMFEqLXTIGEjlu/PqRv1NPdLa1bl2x/GsnLh4QkcBVgCcr2fxTIMYJvpK9swWg0OSZ6zV1dCyfPoO0GB+cviRvX1RW2Z0FePiQkgasAS1C2/6NAjjHhEtlQxmpdFEhIXB4mMuZpYmin7dwZUk0a2bVLuvTS5PqTeWX8PwpkFBMukX3xpZIiRX7DoEBCKjZuDMHr8HAI3IaHw/2sBN5SCLD37Qu30Qh4d/dce1kCb4mrAEtWtv+jQE4RfCMbyhSMUiAhVatWSVu3hhHTrVuzGczm4UNCEvKSKpQZZfo/CuQYwTfSV7ZgtFGBhOg1MzFqZQqy6kMePiR0GlcBlqBs/0eBHCPnG+kbHQ2z8ePBaPyNZGSkWNW6KJDQWWU7n0pgejpMrpycDKkmg4ME3gtw3gOZ0iznm+Ab6SMYRTvVjgAODS28z/mEouH/KJApBN9VBN9ASbDqAwAgRQTfVQTfQImwlCMAICUsNQigXFj1AQCQUQTfAIqFVR8AABn2uLQ7AABt1WgpRym09/ez6gMAIDUE3wCKZWAgLKsWX90hCsD7+0M7AAApIfgGUCxm9Ue2G7UDAJAgcr4BAACAhKQSfJvZajO7zswmqrdHNtjvUTM7UP26Ota+1sy+Vn38XjN7QnK9BwAAAJYnrZHvXZI+7+7rJX2+er+eX7j7hurXq2Ptfy1pqPr4n0ja2tnuAgAAACuXVvB9nqQrqt9fIanlGVBmZpJeLOkzy3k8AAAAkJa0Jlz+urvfK0nufq+ZPbXBfoeb2bikRyS9393HJB0l6afu/kh1n4OSejveYwAosUpF2rtXmpiQ1q+XBgelnp60ewUA+dOx4NvMrpf0G3U2vXsJP+Y4d58ys2dKusHMvi3p53X2a1g1w8wuknSRJB133HFLeGoAgCTt3y9t2hTqFc3MSN3d0o4d0r590saNafcOAPKlY2kn7v5Sd/+tOl+flfTfZnasJFVv72vwM6aqt9+X9EVJp0j6saSnmFn0wWGNpKkm/bjM3fvcve+YY45p2+sDgDKoVELgXamEwFsKt1H79HS6/UNnVSrSnj3Szp3htlJJu0cl4S6Nji6syNuoHbmSVs731ZIurH5/oaTP1u5gZkea2ROr3x8t6UxJt7m7S/qCpAuaPR4AsHJ794YR73pmZ8N2FNP+/VJvr7R9u7R7d7jt7Q3t6LCxMWnzZunii+cCbfdwf/PmsB25lVbw/X5JZ5vZhKSzq/dlZn1mtqe6zwmSxs3sFoVg+/3uflt1205JO8xsUiEH/KOJ9h4ASmJiYm7Eu9bMjDQ5mWx/kAyueKRsYEDatk0aHp4LwC++ONzfto1KvTmXyoRLdz8k6SV12scl/X71+xslndjg8d+X9LxO9hEAECZXdnfXD8C7u6V165LvExpr18TYVq54bM3wIr+5nyBsJg0Nhe+Hh8OXFALvoaGwHbllXqK8ob6+Ph8fH0+7GwCQG5VKSDWol+vb0yNNTUmrViXfLyxUb2JsV9fyJsbu3BlSTRrZtUu69NKV9bdT2nkcUuceOh+ZnSXwzgkzu9nd++pto7w8AKChnp4QtPT0hCBGCrdRO4F3NrQ7TSS64lFPlq94FCpdJko1iYvngCO3CL4BAE1t3BhGuIeHw4jn8HC4n7tRxAJr98TYwcH5A65xXV1hexYVZoJwbY737OzCHHDkVlpFdgAAObJqVbZzfMuu3RNjoysbjdI3snrFI9MThN3DKiUDA/NTR+q1j43NBd5Rjnc8B7y/Xzr//ORfA9qCkW8AAHKuE2kiebzikel0maUsHzgwII2MzJ9cGQXgIyOsdpJzTLgEACABnVyBg4mxQaaPQ20qydDQwvtMpiyMZhMuSTsBAKDD6q3AsWNH+1bgyGuaSLtl+jiwfCCqGPkGAKCDkhyNnZ4Oo+uTkyHFYnCwPIF3XKaPA8sHlgJLDQJAp7hLo6MLVx9o1F4UZX3dkSW8/iRX4Igmxl56abjNTMCZsMweB5YPhAi+AWBlljKJqkjK+rojS3j9mV6BA8lh+UBUkfMNACsxMDD3BiotnERV1FUJyvq6I0t4/dEKHPUC8NRX4EByWD4QVeR8A8BKxUe0ImWYRFXW1x1p8fVnegWOuKWsQ42l4/iWSrOcb4JvAGiHsk6iWuHr7uTye4lo8fXXW+0kWoEjM+tmj46GlJn4B4j4B4yREUZmgRYx4RIAOqmsk6hW+Lr37w8jwtu3S7t3h9ve3tCeC0t4/bkoWBNPpYleR5lSiYCkuHtpvk477TQHgLaanXXfts1dCrf17hfRCl/3z3/u3tMTdq/96ulxr1QSeh3LVdTfe/x1RF95fj1ASiSNe4N4NPWAOMkvgm8AbTcysjBAiQcwIyPp9q9TVvi6/+mf3Lu76wff3d3ue/Yk8BpWosi/99nZ+b8QAm9gyZoF36x2AgArMTAQcmHjk6WiVQz6+4t7qX6Frzv3y+8V9ffuDVJpyjKJFkgAOd8AsBJmYRJabWDSqL0oVvi6o+X36snF8ntF/L1HgTfrUAMdRfANAEjc4OD8RULiurrCdiSs0TrUUQBe9MJJQEIIvgEAievpCcvs9fTMjYB3d8+1Z2Ld67KJUmniKSZRAB6l2ABYMdb5BtA+ThEJLM30dFjne3IypJoMDhJ4A8i/Zut8M+ESQPuMjVGkA0uyapW0dWvavQCA5BB8A2ifeJEOKQTgFOkAAOAxBN8A2ifKD5VCwB0F4fGRcAAASoycbwDt5z5/KYvZWQJvAEBpNMv5ZrUTAO3VqEhHiT7oAwDQCME3gPahSAcAAE2RyskB7wAADWNJREFU8w2gfRoV6ZBCe38/q50AAEqN4BtA+0RFOuLreUcBeH8/q50AAEqP4BtA+5jVH9lu1A4AQMmQ8w0AAAAkhOAbAAAASAjBNwAAAJAQgm8AAAAgIQTfAAAAQEIIvgEAAICEEHwDAAAACSH4BgAAABJC8A0AAAAkhOAbAAAASAjBNwAAAJAQgm8AANB57tLoaLhtpR0oKIJvAADQeWNj0ubN0sUXzwXa7uH+5s1hO1ACqQTfZrbazK4zs4nq7ZF19jnLzA7Evh4ys4Hqto+Z2Q9i2zYk/yoAAEDLBgakbduk4eG5APzii8P9bdvCdqAE0hr53iXp8+6+XtLnq/fncfcvuPsGd98g6cWSHpR0bWyXP4m2u/uBRHoNAHnFJX+kzUwaGpoLwLu65gLvoaGwHSiBtILv8yRdUf3+CkmLfdy9QNLn3P3BjvYKAIqKS/7IgigAjyPwRsmkFXz/urvfK0nV26cusv8WSVfWtL3PzL5lZkNm9sROdBIACoNL/siC6LyLi38gBEqgY8G3mV1vZrfW+TpviT/nWEknSrom1nyJpGdLOl3Sakk7mzz+IjMbN7Px+++/fxmvBAAKoOyX/Em7SV/tB77Z2YUfCIESME/hZDez70l6kbvfWw2uv+juv9lg322SnuvuFzXY/iJJ73L3Vy72vH19fT4+Pr6CngNAzrmHwDsyO1v8wFsKAfbmzfM/bMSDwZER6fzz0+5lsfE7QImY2c3u3ldvW1ppJ1dLurD6/YWSPttk39erJuWkGrDLzEwhX/zWDvQRAIqlzJf8SbtJ38BACLDjV1qiKzIjI/wOUBppjXwfJelTko6T9F+SXuvuD5hZn6Q/dPffr+53vKSvSHq6u8/GHn+DpGMkmaQD1cdML/a8jHwDKK3aYHNoaOH9oo+Ax49BpCyvHUCimo18pxJ8p4XgG0Bpcck/KGvaDYBEZTHtBACQJC75lzvtBkBmEHwDQBmYhZHt2lHeRu1Fw0obADLicWl3AACAjhsbW5jfHhV7GR6W+vvLkXYDIHUE3wCA4ovSbgYGFqbd9PeXI+0GQCYQfAMAii9Kr2m1HQA6hJxvAAAAICEE3wAAAEBCCL4BAACAhBB8AwAAAAkh+AYAAAASQvANAAAAJITgGwAAAEgIwTcAAACQEIJvAAAAICEE3wAAAEBCCL4BAACAhBB8AwAAAAkh+AYAAAASQvANAAAAJITgGwAAAEgIwTcAAACQEIJvAAAAICHm7mn3ITFmdr+ku9PuRwKOlvTjtDuRYxy/5ePYLR/HbmU4fsvHsVsZjt/yFfnYPcPdj6m3oVTBd1mY2bi796Xdj7zi+C0fx275OHYrw/FbPo7dynD8lq+sx460EwAAACAhBN8AAABAQgi+i+mytDuQcxy/5ePYLR/HbmU4fsvHsVsZjt/ylfLYkfMNAAAAJISRbwAAACAhBN85ZWarzew6M5uo3h5ZZ5+zzOxA7OshMxuobvuYmf0gtm1D8q8iHa0cu+p+j8aOz9Wx9rVm9rXq4/ea2ROS6336Wjz3NpjZV83sO2b2LTMbjG0r3blnZuea2ffMbNLMdtXZ/sTquTRZPbeOj227pNr+PTN7WZL9zoIWjt0OM7utep593syeEdtW92+4TFo4fm8xs/tjx+n3Y9surP6dT5jZhcn2PH0tHLuh2HG7w8x+GttW6nPPzC43s/vM7NYG283MPlQ9tt8ys1Nj24p/3rk7Xzn8krRb0q7q97sk/fUi+6+W9ICkJ1fvf0zSBWm/jiwfO0nTDdo/JWlL9ft/kPT2tF9T1o6fpGdJWl/9/mmS7pX0lOr9Up17kg6TdKekZ0p6gqRbJD2nZp8/kvQP1e+3SNpb/f451f2fKGlt9ecclvZrytixOyv2f+3t0bGr3q/7N1yWrxaP31skfbjOY1dL+n719sjq90em/ZqydOxq9v8fki6P3S/7ufdCSadKurXB9k2SPifJJD1f0teq7aU47xj5zq/zJF1R/f4KSQOL7H+BpM+5+4Md7VU+LPXYPcbMTNKLJX1mOY8viEWPn7vf4e4T1e+nJN0nqW6xgRJ4nqRJd/++u/9S0lUKxzAufkw/I+kl1XPtPElXufvD7v4DSZPVn1cWix47d/9C7P/aTZLWJNzHLGvl3GvkZZKuc/cH3P0nkq6TdG6H+plFSz12r5d0ZSI9ywF3/5LCgF8j50n6uAc3SXqKmR2rkpx3BN/59evufq8kVW+fusj+W7TwH8P7qpd7hszsiZ3oZEa1euwON7NxM7spSteRdJSkn7r7I9X7ByX1dra7mbOkc8/MnqcwcnRnrLlM516vpB/G7tc7Zx7bp3pu/UzhXGvlsUW21Ne/VWE0LVLvb7hMWj1+r6n+PX7GzJ6+xMcWVcuvv5rqtFbSDbHmsp97i2l0fEtx3j0u7Q6gMTO7XtJv1Nn07iX+nGMlnSjpmljzJZJ+pBAUXSZpp6T3Lq+n2dOmY3ecu0+Z2TMl3WBm35b08zr7FW7JoDafe5+QdKG7z1abC33u1WF12mrPmUb7tPLYImv59ZvZGyX1SeqPNS/4G3b3O+s9vqBaOX7/JulKd3/YzP5Q4QrMi1t8bJEt5fVvkfQZd3801lb2c28xpf6fR/CdYe7+0kbbzOy/zexYd7+3GuDc1+RHvU7SqLv/Kvaz761++7CZ/bOkd7Wl0xnRjmNXTZeQu3/fzL4o6RRJ/6pweexx1RHKNZKm2v4CUtaO42dmR0j6d0n/b/WyYvSzC33u1XFQ0tNj9+udM9E+B83scZJ+TeGSbSuPLbKWXr+ZvVThg2G/uz8ctTf4Gy5TALTo8XP3Q7G7/yTpr2OPfVHNY7/Y9h5m11L+9rZIeke8gXNvUY2ObynOO9JO8utqSdEs4AslfbbJvgty0apBU5TDPCCp7ozkglr02JnZkVE6hJkdLelMSbd5mBHyBYUc+oaPL7hWjt8TJI0q5PR9umZb2c69b0hab2GVnCcovFHXrn4QP6YXSLqheq5dLWmLhdVQ1kpaL+nrCfU7CxY9dmZ2iqR/lPRqd78v1l73bzixnmdDK8fv2NjdV0v6bvX7aySdUz2OR0o6R/OvnhZdK3+3MrPfVJgY+NVYG+fe4q6W9ObqqifPl/Sz6sBMOc67tGd88rW8L4V80M9Lmqjerq6290naE9vveEn3SOqqefwNkr6tEPj8i6RVab+mLB07SS+oHp9bqrdbY49/pkIANCnp/2/vXkKtquI4jn9/GZZE9BRy4iAIC1Ls5SAsQkQQoZDoQQWhAzGxpg2ksgyKCKJAJDN10KBwZANBMYSkJxKXrglFDTTLBxFRSpblv8FexkG9YrfL8d7r9wMbzn6svfZerMP5s87a+78JuOR839MobL/HgOPAQM8y80Lte3RP9n9DN/K1om17gS5gBLi09aVvW9+6vqfsilbua2D++b6XUdh224FDPf3s/bZ9yO/whbScQ/u9BHzV2mkHcGNP2cWtT34LLDrf9zLa2q6trwRePqXcBd/36Ab8DrTfgf10z2MsBZa2/QFWt7YdBG7vKTvu+50ZLiVJkqQ+cdqJJEmS1CcG35IkSVKfGHxLkiRJfWLwLUmSJPWJwbckSZLUJwbfkjQGJLkmyUBbDib5oWd94n84z+IkZ8peSpKHkuxJciLJzCGOmZBkdZLdSQaTfN7Sa0uSzoEZLiVpDKguE+FMgCQrgSNV9eowTrUY+AI4eIZ9g3SJj9afpfwjdO96n1FVJ5JMBX4dxnX8qydjrCSNewbfkjTGJXmcLr31ROBjYDndP5sb6AL2AGvpktHMBN5L8jswq6r+PHmeqtrTzne26qYAB6rqRCuzr+c6FgCrgAnAoaqa1zL8radL+HUEWFJVu5O8CEymS1p1MMki4BVgNl3SoTeqat3/aBZJGpUMviVpDEtyM7AQuLOq/kqyli4V9nfAtVU1vR13ZVX9kuRJYHlVDQyzyneBnUnuoctw+k5VDbSpLGuAu6pqb5Kr2/GrgM+q6t4k84CNdNlQAW4B7q6qY0mWAYeralZLzf1pkm29wb0kjQcG35I0ts0F7gB2tRHrScD3wFZgWpLXgS3AtpGorKr2JZkGzGnLjiQLgauAHVW1tx33cysyG1jQtm1LsjHJZW3f5qo61j7PA25K8nBbvwK4ATD4ljSuGHxL0tgWYH1VPXPajmQGMB94CrgfWDISFbaAeQuwJclPwH3ATqCGuL6h1o+esn1ZVX0wEtcoSaOVbzuRpLFtO/Bgm1t98q0oU5NMBlJVm4DngFvb8b8Blw+3siS3JZnSPl8ETAf2Ah8Bc06++aRn2smHwKNt21xgf1UdPe3E3Uj9siQXt2OnJZk03OuUpNHKkW9JGsOqajDJ88D2FgwfB5YCfwNvp5uLUsDTrcgGYN2ZHrhM8gDwGt2DkFuT7KqqBadUeR3wVnu9YYBPgDVV9UeSJ4DNrc4f6UbdnwU2JPmS7oHLRUPcypvAVGCgTZ85TDeiLknjSqrO9C+hJEmSpJHmtBNJkiSpTwy+JUmSpD4x+JYkSZL6xOBbkiRJ6hODb0mSJKlPDL4lSZKkPjH4liRJkvrE4FuSJEnqk38AtwS6Gw/4ybQAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 864x576 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "positive = data2[data2['Accepted'].isin([1])]\n",
    "negative = data2[data2['Accepted'].isin([0])]\n",
    "\n",
    "fig, ax = plt.subplots(figsize=(12, 8))\n",
    "ax.scatter(positive['Test 1'],\n",
    "           positive['Test 2'],\n",
    "           s=50,\n",
    "           c='b',\n",
    "           marker='o',\n",
    "           label='Accepted')\n",
    "ax.scatter(negative['Test 1'],\n",
    "           negative['Test 2'],\n",
    "           s=50,\n",
    "           c='r',\n",
    "           marker='x',\n",
    "           label='Rejected')\n",
    "ax.legend()\n",
    "ax.set_xlabel('Test 1 Score')\n",
    "ax.set_ylabel('Test 2 Score')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "这个数据看起来可比前一次的复杂得多。特别地，你会注意到其中没有线性决策界限，来良好的分开两类数据。一个方法是用像逻辑回归这样的线性技术来构造从原始特征的多项式中得到的特征。让我们通过创建一组多项式特征入手吧。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>Accepted</th>\n",
       "      <th>Ones</th>\n",
       "      <th>F10</th>\n",
       "      <th>F20</th>\n",
       "      <th>F21</th>\n",
       "      <th>F30</th>\n",
       "      <th>F31</th>\n",
       "      <th>F32</th>\n",
       "      <th>F40</th>\n",
       "      <th>F41</th>\n",
       "      <th>F42</th>\n",
       "      <th>F43</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>0.051267</td>\n",
       "      <td>0.002628</td>\n",
       "      <td>0.035864</td>\n",
       "      <td>0.000135</td>\n",
       "      <td>0.001839</td>\n",
       "      <td>0.025089</td>\n",
       "      <td>0.000007</td>\n",
       "      <td>0.000094</td>\n",
       "      <td>0.001286</td>\n",
       "      <td>0.017551</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>-0.092742</td>\n",
       "      <td>0.008601</td>\n",
       "      <td>-0.063523</td>\n",
       "      <td>-0.000798</td>\n",
       "      <td>0.005891</td>\n",
       "      <td>-0.043509</td>\n",
       "      <td>0.000074</td>\n",
       "      <td>-0.000546</td>\n",
       "      <td>0.004035</td>\n",
       "      <td>-0.029801</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>-0.213710</td>\n",
       "      <td>0.045672</td>\n",
       "      <td>-0.147941</td>\n",
       "      <td>-0.009761</td>\n",
       "      <td>0.031616</td>\n",
       "      <td>-0.102412</td>\n",
       "      <td>0.002086</td>\n",
       "      <td>-0.006757</td>\n",
       "      <td>0.021886</td>\n",
       "      <td>-0.070895</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>-0.375000</td>\n",
       "      <td>0.140625</td>\n",
       "      <td>-0.188321</td>\n",
       "      <td>-0.052734</td>\n",
       "      <td>0.070620</td>\n",
       "      <td>-0.094573</td>\n",
       "      <td>0.019775</td>\n",
       "      <td>-0.026483</td>\n",
       "      <td>0.035465</td>\n",
       "      <td>-0.047494</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>-0.513250</td>\n",
       "      <td>0.263426</td>\n",
       "      <td>-0.238990</td>\n",
       "      <td>-0.135203</td>\n",
       "      <td>0.122661</td>\n",
       "      <td>-0.111283</td>\n",
       "      <td>0.069393</td>\n",
       "      <td>-0.062956</td>\n",
       "      <td>0.057116</td>\n",
       "      <td>-0.051818</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "   Accepted  Ones       F10       F20       F21       F30       F31       F32  \\\n",
       "0         1     1  0.051267  0.002628  0.035864  0.000135  0.001839  0.025089   \n",
       "1         1     1 -0.092742  0.008601 -0.063523 -0.000798  0.005891 -0.043509   \n",
       "2         1     1 -0.213710  0.045672 -0.147941 -0.009761  0.031616 -0.102412   \n",
       "3         1     1 -0.375000  0.140625 -0.188321 -0.052734  0.070620 -0.094573   \n",
       "4         1     1 -0.513250  0.263426 -0.238990 -0.135203  0.122661 -0.111283   \n",
       "\n",
       "        F40       F41       F42       F43  \n",
       "0  0.000007  0.000094  0.001286  0.017551  \n",
       "1  0.000074 -0.000546  0.004035 -0.029801  \n",
       "2  0.002086 -0.006757  0.021886 -0.070895  \n",
       "3  0.019775 -0.026483  0.035465 -0.047494  \n",
       "4  0.069393 -0.062956  0.057116 -0.051818  "
      ]
     },
     "execution_count": 18,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "degree = 5\n",
    "x1 = data2['Test 1']\n",
    "x2 = data2['Test 2']\n",
    "\n",
    "data2.insert(3, 'Ones', 1)\n",
    "\n",
    "for i in range(1, degree):\n",
    "    for j in range(0, i):\n",
    "        data2['F' + str(i) + str(j)] = np.power(x1, i-j) * np.power(x2, j)\n",
    "\n",
    "data2.drop('Test 1', axis=1, inplace=True)\n",
    "data2.drop('Test 2', axis=1, inplace=True)\n",
    "\n",
    "data2.head()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "现在，我们需要修改第1部分的成本和梯度函数，包括正则化项。首先是成本函数："
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# regularized cost（正则化代价函数）\n",
    "$$J\\left( w \\right)=\\frac{1}{m}\\sum\\limits_{i=1}^{m}{[-{{y}^{(i)}}\\log \\left( {{h}}\\left( {{x}^{(i)}} \\right) \\right)-\\left( 1-{{y}^{(i)}} \\right)\\log \\left( 1-{{h}}\\left( {{x}^{(i)}} \\right) \\right)]}+\\frac{\\lambda }{2m}\\sum\\limits_{j=1}^{n}{w _{j}^{2}}$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [],
   "source": [
    "def costReg(w, X, y, learningRate):\n",
    "    w = np.matrix(w)\n",
    "    X = np.matrix(X)\n",
    "    y = np.matrix(y)\n",
    "    first = np.multiply(-y, np.log(sigmoid(X * w.T)))\n",
    "    second = np.multiply((1 - y), np.log(1 - sigmoid(X * w.T)))\n",
    "    reg = (learningRate /\n",
    "           (2 * len(X))) * np.sum(np.power(w[:, 1:w.shape[1]], 2))\n",
    "    return np.sum(first - second) / len(X) + reg"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "请注意等式中的\"reg\" 项。还注意到另外的一个“学习率”参数。这是一种超参数，用来控制正则化项。现在我们需要添加正则化梯度函数："
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "如果我们要使用梯度下降法令这个代价函数最小化，因为我们未对${{w }_{0}}$ 进行正则化，所以梯度下降算法将分两种情形：\n",
    "\n",
    "\n",
    "\\begin{align}\n",
    "  & 重复\\text{ }直到\\text{ }收敛\\text{ }\\!\\!\\{\\!\\!\\text{ } \\\\ \n",
    " & \\text{     }{{w }_{0}}:={{w }_{0}}-a\\frac{1}{m}\\sum\\limits_{i=1}^{m}{[{{h}}\\left( {{x}^{(i)}} \\right)-{{y}^{(i)}}]x_{_{0}}^{(i)}} \\\\ \n",
    " & \\text{     }{{w }_{j}}:={{w }_{j}}-a\\frac{1}{m}\\sum\\limits_{i=1}^{m}{[{{h}}\\left( {{x}^{(i)}} \\right)-{{y}^{(i)}}]x_{j}^{(i)}}+\\frac{\\lambda }{m}{{w }_{j}} \\\\ \n",
    " & \\text{          }\\!\\!\\}\\!\\!\\text{ } \\\\ \n",
    " & 重复 \\\\ \n",
    "\\end{align}\n",
    "\n",
    "对上面的算法中 j=1,2,...,n 时的更新式子进行调整可得： \n",
    "${{w }_{j}}:={{w }_{j}}(1-a\\frac{\\lambda }{m})-a\\frac{1}{m}\\sum\\limits_{i=1}^{m}{({{h}_{w }}\\left( {{x}^{(i)}} \\right)-{{y}^{(i)}})x_{j}^{(i)}}$\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [],
   "source": [
    "def gradientReg(w, X, y, learningRate):\n",
    "    w = np.matrix(w)\n",
    "    X = np.matrix(X)\n",
    "    y = np.matrix(y)\n",
    "\n",
    "    parameters = int(w.ravel().shape[1])\n",
    "    grad = np.zeros(parameters)\n",
    "\n",
    "    error = sigmoid(X * w.T) - y\n",
    "\n",
    "    for i in range(parameters):\n",
    "        term = np.multiply(error, X[:, i])\n",
    "\n",
    "        if (i == 0):\n",
    "            grad[i] = np.sum(term) / len(X)\n",
    "        else:\n",
    "            grad[i] = (np.sum(term) / len(X)) + (\n",
    "                (learningRate / len(X)) * w[:, i])\n",
    "\n",
    "    return grad"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "就像在第一部分中做的一样，初始化变量。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [],
   "source": [
    "# set X and y (remember from above that we moved the label to column 0)\n",
    "cols = data2.shape[1]\n",
    "X2 = data2.iloc[:,1:cols]\n",
    "y2 = data2.iloc[:,0:1]\n",
    "\n",
    "# convert to numpy arrays and initalize the parameter array w\n",
    "X2 = np.array(X2.values)\n",
    "y2 = np.array(y2.values)\n",
    "w2 = np.zeros(11)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "让我们初始学习率到一个合理值。如果有必要的话（即如果惩罚太强或不够强）,我们可以之后再折腾这个。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {},
   "outputs": [],
   "source": [
    "learningRate = 1"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "现在，让我们尝试调用新的默认为0的$w$的正则化函数，以确保计算工作正常。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.6931471805599454"
      ]
     },
     "execution_count": 23,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "costReg(w2, X2, y2, learningRate)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([0.00847458, 0.01878809, 0.05034464, 0.01150133, 0.01835599,\n",
       "       0.00732393, 0.00819244, 0.03934862, 0.00223924, 0.01286005,\n",
       "       0.00309594])"
      ]
     },
     "execution_count": 24,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "gradientReg(w2, X2, y2, learningRate)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "现在我们可以使用和第一部分相同的优化函数来计算优化后的结果。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(array([ 0.53010248,  0.29075567, -1.60725764, -0.58213819,  0.01781027,\n",
       "        -0.21329508, -0.40024142, -1.37144139,  0.02264304, -0.9503358 ,\n",
       "         0.0344085 ]),\n",
       " 22,\n",
       " 1)"
      ]
     },
     "execution_count": 25,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "result2 = opt.fmin_tnc(func=costReg, x0=w2, fprime=gradientReg, args=(X2, y2, learningRate))\n",
    "result2"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "最后，我们可以使用第1部分中的预测函数来查看我们的方案在训练数据上的准确度。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 26,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "accuracy = 78%\n"
     ]
    }
   ],
   "source": [
    "w_min = np.matrix(result2[0])\n",
    "predictions = predict(w_min, X2)\n",
    "correct = [1 if ((a == 1 and b == 1) or (a == 0 and b == 0)) else 0 for (a, b) in zip(predictions, y2)]\n",
    "accuracy = (sum(map(int, correct)) % len(correct))\n",
    "print ('accuracy = {0}%'.format(accuracy))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "虽然我们实现了这些算法，值得注意的是，我们还可以使用高级Python库像scikit-learn来解决这个问题。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 27,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "LogisticRegression()"
      ]
     },
     "execution_count": 27,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "from sklearn import linear_model#调用sklearn的线性回归包\n",
    "model = linear_model.LogisticRegression(penalty='l2', C=1.0)\n",
    "model.fit(X2, y2.ravel())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 28,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.6610169491525424"
      ]
     },
     "execution_count": 28,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "model.score(X2, y2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "这个准确度和我们刚刚实现的差了好多，不过请记住这个结果可以使用默认参数下计算的结果。我们可能需要做一些参数的调整来获得和我们之前结果相同的精确度。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "## 参考\n",
    "[1] Andrew Ng. Machine Learning[EB/OL]. StanfordUniversity,2014.https://www.coursera.org/course/ml\n",
    "\n",
    "[2] 李航. 统计学习方法[M]. 北京: 清华大学出版社,2019.\n",
    "\n"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.8.8"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}
